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I.  INTRODUCTION 


In  the  November,  1966,  issue  of  Technometrics,  a  paper  'A  Review  of  Response  Surface 
Methodology:  A  Literature  Survey'  by  Hill  and  Hunter  appeared.  That  particular  review  empha¬ 
sized  practical  applications  in  the  chemical  and  processing  fields  and  featured  an  excellent  bibli¬ 
ography.  In  the  December,  1975,  issue  of  Biometrics,  'A  Review  of  Response  Surface  Methodology 
from  a  Biometric  Viewpoint'  by  Mead  and  Pike  appeared.  Emphasis  was  put  on  biological  appli¬ 
cations,  and  a  much  broader  definition  of  RSM  was  made  than  that  of  the  Hill  and  Hunter  paper. 

In  the  Hill  and  Hunter  review,  the  authors  state  that  RSM  began  with  the  work  by  Box  and 
Wilson  (1951).  In  the  Mead  and  Pike  paper,  they  move  back  the  origin  of  RSM  to  include  use  of 
'response  curves'  dating  back  into  the  1930's.  Such  items  as  probit  analysis,  the  use  of  Mitcherlich 
response  equations,  and  the  early  work  in  factorial  arrays  by  Yates  (1935)  are  mentioned  as  prob¬ 
able  motivators  of  the  work  by  Box  and  co-workers.  In  reality  there  are  diverse  notions  as  to  what 
is  the  proper  list  of  topics  to  be  included  under  the  response  surface  umbrella.  This  is  not  surpris¬ 
ing.  The  term  response  surface  analysis  could  be  taken  as  any  analysis  dealing  with  a  fitted  function, 
thus  accommodating  a  large  collection  of  techniques.  Others  view  the  term  as  implying  a  much 
smaller  collection  of  tools.  This  confusion  is  due  in  part  to  the  very  general,  perhaps  even  unfor¬ 
tunate,  name  response  surface.  The  name  implies  something  which  is  broad  and  all  consuming. 
Yet  it  is  interesting  that  all  too  many  subject  matter  scientists,  many  of  which  are  experienced  in 
the  use  of  statistical  tools,  are  not  aware  of  the  term  response  surface  analysis  or  of  the  problems 
it  addresses. 

There  is  no  question  that  the  motivation  for  the  work  by  Box  et  al.  was  the  general  and 
perhaps  ancient  problem  of  planning  and  analyzing  experiments  in  a  search  for  desirable  conditions 
on  a  set  of  controllable,  or  design,  variables,  desirable  often  being  defined  more  strictly  as  those 
conditions  which  give  rise  to  optimum  response.  We  often  refer  to  this  problem  in  a  more  general 
setting  as  one  of  exploring  an  experimental  region,  or  exploring  a  response  surface.  This  gives  rise 
to  the  inevitable  applications  in  the  chemical  and  process  fields  and  the  historical  attraction  to  the 
subject  by  Technometrics  readers.  We  wish  to  make  it  clear  to  the  present  Technometrics  audience 
that  it  is  not  our  intention  to  travel  the  broad  base  of  experimental  design  or  to  review  the  many 
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aspects  of  model  building.  For  example,  we  wish  to  minimize  any  overlap  with  the  recent 
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Technometrics  reviews  by  Steinberg  and  Hunter  (1984)  in  experimental  design  and  Hocking  (1983) 
in  regression  analysis,  though  some  overlap  with  both  will  be  necessary.  ^Ir-intont  is  to  provide 
a  review  of  important  developments  in  response  surface  methodology,  with  the  definition  being 

confined  to  that  of  a  collection  of  tools  in  design  or  data  analysis  that  enhance  the  exploration  of 
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a  region  of  design  variables  in  one  or  more  responses.  -W»  wiH  emphasize  developments  in  the 
statistics  literature  that  have  appeared  since  1966  but  wjll  make  reference  to  some  pre- 1966  work 


where  historical  perspective  is  necessary. 
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In  terms  of  subject  matter  applications,  mt  will  not  at  all  confine  ourselves  to  applications  in 

i 

the  chemical  and  processing  fields.  Wr  feel  that  one  of  the  more  important  items  to  report  is  the 
breadth  of  application  of  RSM  that  has  evolved  in  the  last  20  years.  Applications  have  expanded 
to  areas  such  as  operations  research,  nuclear  energy,  defense  systems  research,  cancer  chemotherapy, 
and  many  others.  We-shall  reviewvthe  spectrum  from  theoretical  developments  to  practical  notions 
of  current  software  availability  in  both  the  design  and  analysis  aspects  of  RSM.  Design  and  analysis 
will  be  treated  separately  though  the  user  must  treat  the  two  together  in  his  total  analysis.  The 
original  intent  by  Box  and  company  was  to  highlight  a  strategy--both  experimental  and  analytic. 
Yet  precious  little  research  employs  a  simultaneous  consideration  of  both  of  these  aspects.  Finally, 
attcmptlu  offer  suggestions  regarding  the  important  areas  for  future  research  and  future  impact 


of  RSM. 
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The  experimental  strategy  and  analysis  in  RSM  revolves  around  the  assumption  that  a  re¬ 
sponse  ii  is  a  function  of  a  set  of  design  variables  x„  jc„  ... ,  x ,  and  that  the  function  can  be  ap¬ 
proximated  in  some  region  of  the  x’s  by  a  polynomial  model.  Prominent  among  the  models 
considered  are  the  first  order  model 


*1  “  Po  +  Ml  +  -  +  M* 


(ID 


or  the  second  order  model 
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Also  in  certain  instances  in  the  sequel  we  will  deal  specifically  with  nonlinear  models.  In  what 
follows,  response  surface  design  is  treated,  with  some  historical  perspective  given,  followed  by  many 
specialized  areas  such  as  mixture  designs,  design  robustness,  alphabetic  optimal  design  theory,  de¬ 
signs  for  nonlinear  models,  and  many  others.  This  is  followed  by  a  treatment  of  response  surface 
analysis  featuring  multiple  responses,  sampling  properties  of  optima,  and  analysis  for  models  other 
than  polynomials.  The  paper  then  features  subject  matter  applications,  with  emphasis  on  industrial 
usage.  Finally,  software  availability  is  reviewed  and  potent  d  areas  for  future  directions  and  re¬ 
search  are  discussed. 

II.  RESPONSE  SURFACE  DESIGN 

In  this  section  we  deal  with  several  aspects  of  response  surface  designs,  beginning  with  the 
'where  we  were'  status  of  response  surface  design  in  the  mid  60's.  We  consider  three  very  general 
areas: 

i)  design  classes 

ii)  optimality  criteria  and  choice  of  specific  designs  that  achieve  certain  important  prop¬ 
erties 

iii)  designs  which  accomplish  special  design  goals  or  have  special  design  features. 

Item  (i)  may  be  viewed  as  design  families,  e.g.,  composite  designs,  classes  of  three-level  designs,  etc. 
Item  (ii)  reviews  criteria  from  which  choices  of  design  parameters  in  (i)  can  be  made;  for  example, 
the  choice  of  a  composite  design  that  is  rotatable,  or  the  choice  of  a  three-level  design  that  is  D- 
optimal.  In  item  (iii)  we  consider  less  general  topics  such  as  designs  for  nonlinear  models,  robust 
designs,  designs  for  mixture  problems,  designs  for  estimating  slope,  etc.  One  may  historically  view 
the  development  of  RS\1  designs  as  taking  on  these  three  general  areas  of  concentration.  The  user 
has  relied  heavily  on  the  application  of  properties  «jtd  criteria  in  (ii)  and  special  features  in  (iii)  to 
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design  families  in  (i)  that  can  be  used  in  practice.  In  these  areas  much  progress  has  been  made  since 
the  mid  60's. 

21  Status  of  Response  Surface  Design  in  Mid  1960's 

Among  the  many  important  works  that  preceded  the  Hill  and  Hunter  review,  three  are  note* 
worthy  here.  The  first  was  the  benchmark  paper  by  Box  and  Wilson  (1951)  in  which  the  notion 
of  composite  designs  was  introduced.  The  introduction  of  the  'axial  portion*  to  augment  a  two- 
level  factorial  array  was  done  to  allow  for  efficient  estimation  of  quadratic  terms  in  the  second  order 
model  of  equation  (1.2).  This  class  of  designs  allows  flexibility  in  the  choice  of  the  axial  parameter 
and  number  of  center  runs.  This  would  later  prove  very  beneficial  as  more  sophisticated  properties 
and  criteria  evolved.  More  details  regarding  the  composite  designs  are  given  in  section  2.2.  Sharing 
the  sphere  of  influence  with  Box- Wilson  is  the  Box-Hunter  paper  (1957)  in  which  the  notion  of 
rotatability  was  introduced.  The  property  of  mutability  requires  that  the  variance  of  a  predicted 
value  remain  constant  at  points  that  are  equidistant  from  the  design  center.  Rotaubility  was  and 
remains  an  important  design  property.  The  reader  can  gain  insight  into  its  importance  by  observing 
that  in  much  of  the  RSM  design  research  that  appeared  in  the  late  60  s  and  70's,  investigators  began 
with  the  presumption  of  mutability.  No  one  doubts  the  elegance  of  the  work  that  led  to  the 
rotaubility  property.  Most  prudent  consultants  in  RSM  do  not  recommend  designs  that  have 
considerable  deviation  from  rotaubility.  The  influence  of  the  Box-Hunter  work  was  enhanced 
greatly  by  the  fact  that  the  conditions  for  rotatability  are  so  easily  attained  in  the  case  of  first  and 
second  order  designs. 

A  third  very  influential  piece  of  work  of  the  50's  and  60  s  appeared  in  Box  and  Draper  (1959, 
1963).  It  was  an  introduction  of  the  notion  of  robustness  of  an  RSM  design  to  model  misspecifi- 
cation,  though  the  word  'robustness'  was  not  used  in  this  regard  until  later.  A  mean  squared  error 
type  of  design  criterion  was  introduced  which  accounts  for  bias  introduced  when  the  fitted  RSM 
model  is  a  polynomial  of  Iowct  degree  than  the  'true'  model.  Convincing  arguments  are  made  that 
not  only  can  bias  not  be  ignored  but  that  if  there  is  even  a  modest  amount  of  misspecification,  the 
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approach  of  choosing  the  design  that  'protects*  against  bias  docs  not  deviate  substantially  from  that 
of  minimization  of  mean  squared  error. 

The  fundamental  philosophy  of  the  Box* Draper  work  centered  around  the  consideration  of 
the  average  mean  squared  error 

Jm~ f  fa 

where  y(^)  is  the  fitted  polynomial  of  order  </,  and  g{g)  is  a  model  of  order  di>  dx  containing  un¬ 
known  parameters  and  is  regarded  as  the  'true'  mean  response  or  at  least  it  can  be  viewed  as  the 
model  that  one  chooses  to  protect  against.  Here,  R  is  the  region  of  interest  in  the  design  variables, 
K  is  the  reciprocal  of  the  volume  of  R,  N  is  the  total  number  of  observations,  and  a2  is  the  error 
variance.  The  expression  J  divides  into  the  sum  of  the  variance  and  squared  bias  averaged  over  the 
region  R.  The  notion  of  a  minimum  bias  design  was  introduced.  In  the  formulation  presented  here, 
the  minimum  J  design  cannot  be  achieved  and  evidence  was  put  forth  that  a  strategy  of  design 
choice  that  minimizes  the  bias  portion  of  J  is  an  effective  approach  across  a  broad  range  of  model 
misspecification.  The  importance  of  the  Box-Draper  work  lies  in  the  fact  that  much  of  the  more 
specialized  RSM  design  work  that  followed  into  the  late  60's  and  70's  was  flavored  by  the  now  es¬ 
tablished  need  to  consider  model  underspecification  in  any  serious  attempt  in  developing  optimal 
designs. 

While  much  of  the  foundation  of  very  fundamental  work  was  established  by  that  cited  in  the 
foregoing,  the  appetite  of  the  practical  user  of  RSM  had  been  satisfied  by  the  establishment  of 
families  of  useful  experimental  designs  for  first  order  and  second  order  models.  In  the  first  order 
case,  the  need  for  orthogonal  designs  was  rr\otivated  in  the  Box-Wilson  paper,  Box  (1951),  and  an 
excellent  text  Design  and  Analysis  of  Industrial  Experiments  by  Davies  (1954).  Specific  design 
classes,  two-level  factorial  and  fractional  factorial  designs  had  been  discussed  at  length  in  Box  and 
Hunter  (1961a,  1961b).  Simplex  designs  and  Plackctt-Burman  (1946)  designs  were  available  as 
economical  first  order  orthogonal  designs.  For  second  order  models,  many  subject  matter  scientists 
and  engineers  have  a  working  knowledge  of  the  family  of  central  composite  designs  and  a  class  of 
special  three-level  designs  by  Box  and  Behnken  (1960).  Another  important  pre-1966  contribution 


came  from  Manley  (1959),  in  which  an  effort  was  made  to  reduce  the  number  of  runs  in  the  com* 
posite  design  from  those  that  had  earlier  been  suggested.  It  also  allowed  for  more  flexibility  in 
choosing  two*levcl  fractions  for  the  cube  portion  of  the  composite  design. 

The  majority  of  response  surface  real-life  applications  involve  the  use  of  a  fust  or  second  order 
model.  Where  a  second  order  model  is  inappropriate,  a  nonlinear  model  is  generally  used.  How¬ 
ever,  there  are  instances  in  which  a  third  order  model  is  successful,  and  thus  there  is  a  need  for  third 
order  designs.  By  the  mid  I960's,  a  majority  of  the  work  dealing  in  third  order  designs  essentially 
involved  searching  for  designs  that  possess  the  property  of  rotatability.  Included  are  papers  by 
Draper  (1960,  1961,  1962),  Gardiner  et  al.  (1959),  and  Herzberg  (1964). 

In  the  following  sections  we  provide  reviews  of  what  we  feel  are  important  RSM  design  top¬ 
ics,  with  the  chronology,  for  the  most  part,  beginning  in  the  mid  I960's.  We  begin  with  specific 
families  of  second  order  designs,  with  some  emphasis  put  on  alternatives  to  the  central  composite 
designs,  though  the  latter  remains  the  most  flexible  and  thus  most  popular  family  of  second  order 
designs  used  in  practice. 

2.2  The  Central  Composite  Design  and  Alternate  Design  Plans 

The  central  composite  design  (c.c.d.)  is  a  special  case  of  the  class  of  composite  designs  in 
which  the  construction  consists  of  three  portions: 

i)  The  2*  vertices  of  a  cube  (or  a  fraction  of  these  vertices) 

ii)  The  2k  vertices  of  a  cross- polytope  or  'star',  or  axial  portion  with  parameter  a 

iii)  A  chosen  number,  %  of  center  runs. 

The  points  in  (i)  form  a  two-level  factorial  or  fractional  factorial  anay.  The  points  in  (ii)  essentially 
form  an  augmentation  which  represents  a  one-factor-al-a-time  portion  designed  to  provide  esti¬ 
mation  of  the  {5,„  the  pure  quadratic  terms  in  equation  (1.2).  Center  runs  in  portion  (iii)  obviously 
provide  an  internal  estimate  of  pure  error  variance.  The  choice  of  the  number  of  center  runs  is  very 
important  to  the  practitioner.  Figure  1  provides  an  illustration  of  the  c.c.d.  for  k  =  3  design  vari¬ 
ables. 
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Figure  1.  A  Central  Composite  Design  for  k  *  3  Design  Variables 

The  design  levels  are  centered  and  scaled  to  design  units  with  ±  1  being  the  levels  in  the  factorial 
portion  and  the  parameter  a  being  the  axial  distance  as  indicated  in  Figure  1.  In  large  part,  the 
appeal  of  the  c.c.d.  is  derived  from  its  flexibility  and  due  to  the  fact  that  it  is  a  natural  design  plan 
for  sequential  experimentation,  since  the  factorial  portion  and  center  runs  served  as  a  preliminary 
phase  from  which  one  may  fit  a  first  order  model  and  yet  gain  evidence  regarding  the  importance 
of  pure  quadratic  contribution.  If  pure  quadratic  terms  are  found  to  be  significant  on  the  basis  of 
the  preliminary  analysis,  the  design  may  be  augmented  with  axial  runs  for  the  completion  of  the 
second  order  fit.  For  a  fundamental  look  at  the  central  composite  design,  the  reader  is  referred  to 
Davies  (1954),  Myers  (1976),  Box  et  al.  (1978),  and  Cochran  and  Cox  (1957). 

Much  of  the  flexibility  of  the  central  composite  design  is  derived  from  the  user's  choice  of  a 
and  n 0.  There  are  many  criteria  used  for  the  choice  of  these  parameters. 

Choice  of  a 

For  a  >  k'  \ the  axial  points  are  further  from  the  design  center  than  the  factorial  points.  A 
value  a  =  (F)1'*,  where  F  is  the  number  of  factorial  points,  guarantees  rotatability.  Rotatability  in 
the  c.c.d.  was  prominent  in  the  early  response  surface  design  research  in  the  late  50's  and  60's. 
However,  other  criteria  such  as  robustness  to  model  misspecificahon,  quality  estimation  of  the 
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slope,  robustness  to  errors  in  design  levels,  robustness  to  outliers,  generalized  variance  of  model 
coefficients,  and  orthogonal  blocking  are  among  the  important  considerations  for  choice  of  a  con- 
sidered  by  researchers.  All  of  these  are  discussed  in  later  sections  where  these  topics  or  criteria  are 
introduced.  In  many  eases,  of  course,  the  property  chosen  by  the  experimenter  cannot  be  satisfied 
with  the  rotatable  value  a  ■  (f)'1*.  It  is  generally  considered  that  while  exact  rotatability  is  not 
necessary  in  practice,  a  substantial  deviation  is  not  wise. 


Choke  of  /Tq 


The  choice  of  the  number  of  center  runs  is  a  vital  ingredient  in  the  successful  use  of  the  central 
composite  design  for  the  construction  of  a  second  order  response  surface.  Centf  runs  have  long 
been  considered  important  as  a  source  of  replication  error  in  the  analysis.  From  Box  and  Hunter 
(1957),  it  becomes  apparent  that  the  early  research  into  the  proper  choice  of  r \  presumed  the  design 
in  question  is  indeed  rotatable.  The  recommendation  was  to  choose  ^  for  which  uniform  infor¬ 
mation,  or  uniform  precision  is  achieved;  this  property  provides  approximately  a  uniform  value  of 
Var  y(jt)  inside  a  sphere  of  a  specific  radius.  The  purpose  of  this,  of  course,  is  to  produce  a  stability 
in  the  prediction  variance  in  the  region  of  interest  of  the  design  variables  or,  rather,  to  control  the 
prediction  variance  in  such  a  way  that  equal  precision  in  prediction  is  achieved  everywhere  in  some 
portion  of  the  design  region.  It  is  a  natural  extension  to  the  rotatability  property,  the  former  re¬ 
quiring  Var  yU)  to  be  constant  on  a  sphere.  However,  there  has  long  been  confusion  among  usen 
about  what  constitutes  this  subregion.  See  Draper  (1982).  The  original  Box- Hunter  definition  of 

a 

uniform  precision  requires  n^  to  be  chosen  such  that  Var  y{s)  is  approximately  constant  inside  a 

sphere  of  radius  1.0,  with  the  metric  of  design  units  requiring  scaling  so  that  the  second  moment 

of  the  design  variables  is  unity.  In  the  metric  of  the  standard  ±  1  scaling  for  the  factorial  portion 

of  the  c.c.d.,  this  is  equivalent  to  a  radius  p  -  f-=r — tt*- -  •  As  an  example,  if  a  uniform 

-J  r  +  2k  +  r\, 

precision  rotatable  design  for  k  ■  2  is  constructed,  the  axial  parameter  a  **  \fl  and  ^  3  5.  (See 
Myers  (1976)).  The  implication  is  that  the  prediction  variance  is  roughly  constant  inside  a  sphere 
of  radius  p  "  “  0.83  .  This  does  not  seem  particularly  pleasing  in  light  of  the  fact  that,  in 

this  case,  all  design  points  are  at  a  distance  p  =  y/T  from  the  design  center,  and  no  control  is  ex¬ 
erted  oq  Var  y{&)  at  a  distance  which  approaches  the  perimeter  of  the  design  region.  Thus  the 
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property  of  uniform  precision  does  not  guarantee  stability  in  Vor  yOs.)  everywhere  in  the  design 
region,  and  yet  this  fact  is  either  overlooked  by  or  unknown  to  subject  matter  users. 

As  we  have  described  here,  many  of  the  recommendations  regarding  the  number  of  center 
runs  for  a  c.c.d.  have  been  confined  to  the  rotatable  c.c.d.  Lucas  (1977)  offers  alternative  ideas  on 
number  of  center  runs  for  both  the  family  of  composite  designs  and  other  designs.  An  excellent 
general  discussion  of  center  runs  in  response  surface  designs  appears  in  Draper  (1982). 

2.3  Orthogonal  Blocking 

The  central  composite  design  is  the  second  order  family  that  is  rich  in  candidate  designs  that 
block  orthogonally.  The  latter  describes  a  condition  in  which  regression  coefficients  are  orthogonal 
to  block  effects  and  thus  the  analysis  can  be  conducted  with  minimal  effort  free  from  blocks.  Box 
and  Hunter  (1957,  1961a,  1961b)  developed  the  general  conditions  that  give  rise  to  orthogonal 
blocking  in  the  second  order  case.  The  Box-Behnken  design  also  gives  rise  to  orthogonal  blocking 
in  some  situations. 

2.4  Other  Second  Order  Designs 

Though  the  central  composite  remains  the  second  order  design  family  that  is  most  often  used, 
there  are  other  design  classes  that  are  attractive  in  various  circumstances.  We  earlier  made  reference 
to  the  class  of  Box-Behnken  designs.  The  Box-Behnken  designs  can  be  very  useful  in  cases  when 
it  is  important  that  three  levels  be  used.  In  what  follows,  we  review  other  second  order  designs  that 
appeared  in  the  literature  after  the  1960's. 

Hoke  Design 

Hoke  (1974)  developed  a  class  of  economical  second  order  designs  based  on  irregular  fractions 
of  the  3"  factorial  plan.  These  fractions  are  based  on  sets  of  partially  balanced  arrays  and  are 
compared  to  Box-Behnken  and  other  competing  designs. 


Hybrid  Design 


Roquemore  (1976)  developed  a  family  of  second  order  designs  for  k  ■  3,  4,  6,  and  7  that  are 
either  saturated  or  near  saturated  and  have  some  similarities  to  the  centra]  composite  design.  The 
design  for  k  variables  is  constructed  by  augmenting  a  (k  -  1)  dimensional  central  composite  design 
with  an  additional  column  in  the  design  matrix.  The  values  for  the  additional  column  are  chosen 
to  achieve  design  moment  characteristics  that  allow  a  degree  of  orthogonality  that  is  similar  to  that 
of  the  c.c.d.  It  is  of  interest  to  show  an  example.  For  k  ■  3,  we  have  the  following  design  matrix 


*1 

—  1 

-1 

0.6386 

—  1 

1 

O.o386 

1 

-1 

0.6386 

1 

1 

0.6386 

1.1736 

0 

-0.9273 

-1.1736 

0 

-0.9273 

0 

-1.1736 

-0.9273 

0 

1.1736 

-0.9273 

0 

0 

1.2906 

0 

0 

-0.1360 

Note  the  central  composite  type  structure  in  x ,  and  x,  while  constant  levels  in  Xj  reside  at  the  fac¬ 
torial  portion  and  the  axial  portion  of  x,  and  x,.  In  addition,  two  points  enter  along  the  axis  of 
X).  Center  runs  are  also  suggested.  The  constants  in  the  design  are  determined  so  that  all  pure 
second  moments  are  equal  and  all  odd  moments  are  zero. 

Other  Second  Order  Economical  Designs 

The  work  by  Hartley  (1959)  cited  earlier  and  additional  work  by  Westlake  (1965)  and  Draper 
(1985)  produced  a  class  of  designs  called  small  composite  designs.  These  are  central  composite  de¬ 
signs  but  they  are  worthy  of  special  note  because,  like  the  class  of  hybrid  designs,  they  find  con¬ 
siderable  use  when  experimentation  is  costly.  The  designs  allow  considerable  reduction  of  cost 
below  what  is  required  by  the  original  and  perhaps  standard  concept  of  the  c.c.d.,  the  latter  re¬ 
quiring  the  design  to  contain  a  full  2*  or  a  fraction  of  resolution  V  or  more,  combined  with  the  axial 
portion.  The  purpose  of  this  restriction,  of  course,  is  to  allow  orthogonality  among  the  linear  co¬ 
efficients  and  two  factor  interactions  in  the  second  order  model.  The  notion  of  'small'  composite 
results  from  the  fact  that  additional  information  on  linear  coefficients  is  obtained  from  the  axial 
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portion  of  the  design  though  the  interaction  coefficients  are  estimated  only  from  the  factorial  por¬ 
tion.  One  may  take  advantage  of  this  by  using  as  the  factorial  portion  a  fraction  in  which  linear 
coefficients  are  aliased  with  two  factor  interactions.  An  example  of  a  small  composite  design  for 
three  design  variables  is  given  by 


1  -1 
-1  1 
-o  0 
a  0 
0  -a 
0  a 


1  1 
1  1 
0  0 
0  0 
a  0 
a  0 
0  -a 


The  above  defines  a  saturated  second  order  design  in  which  all  coefficients  in  the  standard  second 
order  model  are  estimable.  Draper  (1985)  demonstrated  that  in  some  cases,  Plackett-Bunnan  de¬ 
signs  can  be  used  in  the  factorial  portion. 

The  classes  of  small  composite  designs  and  hybrid  designs  are  by  no  means  the  only  choices 
that  are  available  when  one  seeks  a  saturated  or  near  saturated  second  order  design.  For  example, 
Notz  (1982)  suggested  a  method  of  constructing  designs  that  are  very  efficient  in  terms  of  general¬ 
ized  variance,  i.e.,  D-efficiency.  (For  more  complete  details  on  the  use  of  D-efficiency  and  D- 
optimality,  see  Section  2.S.)  His  designs  are  saturated  and  constructed  from  the  3*  lattice.  For 
saturation,  we  have  S  *  ~(k  +  l)(/c  +  2)  ■  p,  the  number  of  parameters  in  a  second  order  model, 
and  the  design  matrix  is  given  by 


where  D,  is  a  (p  -  k)  *  k  matrix  chosen  as  a  specific  irregular  fraction  of  a  2*,  whereas  is  a 
k  *  k  taken  from  (0, 1)*  with  at  least  one  coordinate  being  zero.  For  example,  the  k  =  3  design 
matrix  is  given  by 


II 


The  designs  are  shown  to  compare  favorably  (in  the  sense  of  D-cfficiency)  with  other  economical 
second  order  designs. 

Box  and  Draper  (1971,  1974)  produce  practical  designs  that  are  saturated,  efficient  from  a 
generalized  variance  standpoint,  and  contain  four  levels  for  each  variable.  An  example  for  k  »  3 
is  as  follows: 


*1 

-1 

-1 

-1 

1 

-1 

-1 

-1 

1 

-1 

-1 

-1 

1 

X 

X 

-1 

X 

-1 

X 

-I 

X 

X 

1* 

1 

1 

1 

n 

1 

1 

1 

y 

where  X  ”  0.1925  and  \i  “  -0.2912.  The  X  and  p  are  chosen  to  minimize  the  generalized  variance 
of  the  regression  coefficients  in  the  second  order  model. 

During  the  70's  and  early  80's,  much  of  the  work  that  resulted  in  the  development  of  new 
second  order  designs  revolved  around  the  concept  of  D-optimality  and  D-effiriency.  Some  of  this 
work  produced  designs  that  are  very  practical  and,  from  the  foregoing,  one  can  see  that  many  useful 
saturated  or  near  saturated  designs  came  from  that  search.  The  consideration  of  D-eflicicncy  and 
D-optimality  also  produced  other  interesting  results  that  allowed  for  a  better  understanding  of  the 
comparison  of  competing  designs.  In  a  very  interesting  and  timely  paper  by  Lucas  (1974),  he  em¬ 
phasizes  that  in  terms  of  generalized  variance,  if  there  is  an  experimental  region  fixed  by  practical 
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considerations  (say  a  cube  or  hypcrsphcrc),  the  optimum  central  composite  design  is  one  in  which 
the  experimental  points  are  moved  to  the  edge  of  the  experimental  region.  If  the  permissible  ex¬ 
perimental  region  is  a  hypercube,  the  factorial  points  should  be  at  the  vertices  of  the  hypcrcubc 
(  ±  1)  and  the  axial  parameter  a  ■  1.0  .  If  the  region  is  a  hypersphere  with  radius  V k  ,  then  the 
optimum  a  ■  s/k. 

2.5  Application  of  Criteria  Dealing  with  Optimal  Design  Theory 

The  important  development  of  optimal  design  theory  in  the  field  of  experimental  design 
emerged  following  World  War  II.  It  was  motivated  by  various  authors  including  Elfving  (1952, 
1955,  1959),  Chemoff  (1953),  Kiefer  (1958,  1959,  1960,  1961,  1962a,  1962b)  and  Kiefer  and 
Wolfowitz  (1959,  1960).  Kiefer,  in  particular,  was  instrumental  in  providing  the  mathematical 
groundwork  for  optimal  design  theory.  Actually,  the  origin  of  this  theory  can  be  traced  back  to  a 
paper  by  Smith  (1918).  Some  early  work  in  this  area  was  done  by  Wald  (1943)  and  Mood  (1946). 

Optimal  design  theory  has  clearly  become  an  important  component  in  the  general  develop¬ 
ment  of  experimental  design  for  the  case  of  regression  models.  It  also  had  an  impact  on  how  re¬ 
searchers  viewed  and  proposed  response  surface  design  criteria  in  the  1970's  and  early  1980's. 
Methods  were  proposed  for  choosing  RSM  design  parameters  and  developing  second  order  eco¬ 
nomical  designs  through  the  use  of  this  theory.  Though  no  one  can  reasonably  deny  contribution 
and  impact  of  optimal  design  theory,  there  are  many  who  feel  as  if  its  use  in  RSM  design  work 
should  be  done  with  extreme  caution.  In  fact,  the  contrast  in  the  views  concerning  this  issue  has 
resulted  in  what  some  refer  to  as  the  'Kiefer  approach'  and  others  refer  to  as  the  "Box  approach'. 
In  Kiefer's  approach  a  design  is  a  probability  measure,  s,  defined  on  a  closed  and  bounded  subset 
X  of  a  Euclidean  space  of  dimension  equal  to  the  number  of  design  variables  in  the  fitted  model. 
Hence,  this  measure  must  satisfy  the  conditions 

r(s)  i  0  for  all  a  e  x 
Ms)  =  l  • 
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In  particular,  a  collection  of  N  points  in  x.  not  necessarily  distinct,  form  a  design  measure.  In  this 
ease 

f  rtf/N  if  ^  is  a  design  point  (i  »  1,  2,  ... ,  m) 

*<*)■  <  .  .  .  (21) 
(  0  otherwise  , 

where  n,  denotes  the  number  of  replications  at  the  /  th  design  point,  and  m  is  the  number  of  distinct 
design  points.  Such  a  design  measure  is  said  to  be  discrete  and  is  denoted  by  .  Discrete  design 
measures  are  the  traditional  designs  considered  in  Box's  approach.  Design  measures  other  than 
discrete  are  said  to  be  continuous.  These  include  design  measures  of  the  form 


*%)  * 


,  X,  if  £/  is  a  design  point  (/  “  1,  2,  ... ,  l) 
\  0  otherwise  , 


(2.2) 


where  )i,  10  and  x  X  ,*  1  ,  and  at  least  one  X,  is  an  irrational  number,  that  is,  it  is  not  expressible 

/-i 

as  a  fraction.  A  continuous  design  measure  is  sometimes  referred  to  as  an  approximate  design 
measure  since  it  is  not  realizable  in  practice,  but  can  be  approximated  fairly  closely  by  a  discrete 
design  of  the  form  given  in  (2.1).  The  latter  design  is,  therefore,  referred  to  as  an  exact  design 
measure. 

Suppose  that  the  fitted  model  is  of  the  form 


y  ■  f(s)0  +  e. 


(2.3) 


where  y  is  the  measured  response  at  a  point  s  e  x,  f  U)  is  a  vector  of  known  functions  of  £,  fi  is 
an  unknown  vector  of  p  parameters,  and  c  is  a  random  error.  If  the  design  measure  is  discrete,  then 
model  (2.3)  can  be  written  in  vector  form  as 

l  -  +  i,  (2  4) 


where  X  is  an  /V  *  p  matrix  whose  /  th  row  consists  of  the  elements  of  f(&)  (/  =  1,  2,  ... ,  N). 
In  general,  for  any  design  measure  s  on  x.  the  information  matrix  of  s,  denoted  by  A f(s),  is  defined 
as  the  p  *  p  matrix  |my  (r)]  ,  where 


ttty(j)  -  Jx  f/(i)f/*)cij(*)  ,  i,j  ■  1,  2 . p  , 


(2.5) 


where  ft(jc)  is  the  /th  element  of  f(£)  (/  ■  1,  2,... ,  p)  .  In  particular,  if  s  is  a  discrete  design 
measure,  DN,  then  M{DN)  “  X’X/N,  where  A"  is  the  matrix  in  (2.4).  In  this  case  if  we  assume  that 
£( t)  *  Q  and  Varft)  »  o1/*,  then  the  variance  of  the  predicted  response  at  a  point  £  is 

Var  J(*)  -  o2f'(a)(r^)",f(2)  .  (2.6) 

Let  v(i,  DN)  denote  the  standardized  prediction  variance  ( N/a })  Var  y{&),  then 

ivrte)(r^)"'f(i)  .  (2.7) 

If  £  in  (2.4)  is  normally  distributed,  then  a  confidence  region  for  £  of  a  given  confidence  coefficient 
has  the  form 


(fi:  (fi  -  G)' AT AT(fi  -  fi)  s  constant}  , 

which  is  an  ellipsoid  centered  at  £,  the  least  squares  estimator  of  £  .  The  volume  of  this  ellipsoid 
is  proportional  to  \X’X\  ~m  ,  where  lATAfl  denotes  the  determinant  of  X"X.  This  volume  can, 
therefore,  be  minimized  if  the  discrete  design  DN  is  chosen  such  that  \X"X\  is  maximum  over  the 
class  of  jV-point  discrete  designs.  A  design  having  this  property  is  said  to  be  a  discrete  D-optimal 
design. 

In  general,  for  any  design  measure  s  defined  on  x>  the  standardized  prediction  variance  de¬ 
scribed  in  (2.7)  will  have  the  form 


vOc,  s)  -  f,(x)A/"1(j)f(a)  (2.8) 

provided  that  the  information  matrix  M{j)  is  nonsingular.  If  H  denotes  the  class  of  all  design 
measures  on  x,  then  a  design  measure  s*  is  D-optimal  if  it  maximizes  I  ,Vf(j)  |  over  H,  that  is, 

l/V/(j*)|  »  sup  I  \f(s)  I  .  (2.9) 

ieH 

IS 
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This  is  an  extension  of  the  concept  of  D-optimality  in  the  discrete  case. 

The  D-optimality  criterion  has  received  a  great  deal  of  attention  as  evidenced  by  the  numer¬ 
ous  articles  written  about  it.  The  review  articles  by  St.  John  and  Draper  (1975),  Ash  and  Hedayat 
(1978),  and  more  recently  Atkinson  (1982)  contain  many  references  on  D-optimality  and  other 
optimal  design  criteria. 

The  objective  function  in  the  D-optimality  criterion,  namely  I  M(s)  I ,  can  be  expressed  in 
terms  of  the  eigenvalues  X,  (r)  of  the  matrix  M(s)  since  I M(s)  I  *  riX,  (s)  ,  where  p  is  the  number 
of  parameters  in  the  model.  Other  functions  of  these  eigenvalues  were  considered  in  the  develop¬ 
ment  of  mainly  continuous  design  measures.  These  include  tr|A/(r)J  **  t\,  (j)  and  X(1)(r),  where 
X(l)(r)  is  the  smallest  eigenvalue  of  M(s).  The  design  criteria  based  on  these  functions  are  called 

I 

A-optimality  and  E-optimality,  respectively.  More  specifically, 

i)  a  design  measure  s*  is  A-optimal  if  it  maximizes  iX^s)  over  H,  the  class  of  all  design 
measures  on  x- 

ii)  a  design  measure  s*  is  E-optimal  if  it  maximizes  X(1)(r)  over  H. 

Another  design  criterion  which  is  closely  related  to  D-optimality  is  the  G-optimality  criterion.  A 
design  measure  s*  is  G-optimal  if  it  minimizes  over  H  the  maximum  variance  function  given  in 
(2.8),  where  the  latter  maximization  is  carried  out  over  the  experimental  region  x.  that  is, 

sup  Kx.-t*)  *  inf  /sup  v(i, s)\,  (2.10) 

aex  \zex  / 

Kiefer  and  Wolfowitz  (1960)  proved  that  the  D-optimality  and  the  G-optimality,  as  defmed  in  (2.9) 
and  (2.10),  respectively,  are  equivalent.  Furthermore,  they  showed  that  these  two  extremum 
problems  are  equivalent  to  yet  another  extremum  problem,  namely 

sup  v(a,  j*)  *  p  .  (2.11) 

sex 


It  is  to  be  noted  that,  in  general, 
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sup  v(jc,  j*)  i  p 

t  n 


(2-12) 


(see  Fedorov  1972,  p.  69).  Hence,  the  maximum  standardized  variance  of  the  predicted  response 
achieves  its  minimum  value  of  p  whenever  the  design  measure  s*  is  D*  or  G-optimal. 

The  above  results  can  be  formalized  in  the  following  theorem,  known  as  the  Equivalence 
Theorem  of  Kiefer  and  Wolfowitz: 

Theorem:  If  H  is  the  class  of  all  design  measures  defined  on  a  subset  %  of  the  ^-dimensional 
Euclidean  space  ( k  is  the  number  of  design  variables  in  the  model),  and  if  p  is  the  number  of  pa¬ 
rameters  in  the  model,  then  the  following  assertions  are  equivalent: 

i)  s*  is  D-optimal 

ii)  s*  is  G-optimal 

iii)  sup  vte.j*)  ■  p. 

Assertion  (iii)  of  the  Equivalence  Theorem  is  very  useful.  It  provides  an  easy  check  on 
whether  a  design  measure  is  D-optimal  by  comparing  the  maximum  of  v{s,  s)  over  %  w>th  the 
number  of  parameters  in  the  corresponding  model.  Furthermore,  it  can  be  effectively  used  to 
construct  a  D-optimal  design  from  a  nondegenerate  discrete  design  4*  (nondegenerate  means  that 
the  corresponding  XX  matrix  is  nonsingular)  by  augmenting  it  sequentially  with  design  points 
chosen  appropriately  (see  Wynn  1970,  Silvey  1980,  Chapter  4).  The  first  point,  *,  is  chosen  such 
that 


,DNq)  =  sup  (2.13) 

in 

By  augmenting  D „0  with  z,  we  get  ,  which,  in  turn,  is  augmented  with  the  point  z}  which  sat¬ 
isfies 


><*2’  DNl)  -  sup  vU,  Dn  ). 
in 


(2.14) 
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This  process  continues  and  results  in  the  sequence  of  designs  DNt)  c  D*,  c  ....  which,  according 
to  the  Equivalence  Theorem,  must  converge  to  a  D-optimal  design  measure  (see  Fedorov  1972, 
Th.  2.5.3;  Wynn  1970,  Th.  1).  Note  that  since  by  (2.12), 

KSj.  0/vM)  ~  P  *  0  i  *“  1.2,...,  (2.15) 

the  sequential  generation  of  a  D-optimal  design  can  be  terminated  whenever  the  difference  in  (2.15) 
is  less  than  some  value  5  >  0  chosen  a  priori.  Thus,  it  is  possible  to  approximate,  with  a  given 
precision,  a  D-optimal  design  measure  with  a  discrete  design.  This  fact  has  in  effect  brought  D* 
optimal  designs  out  of  the  realm  of  theory  into  factual  existence.  Nalimov  et  al.  (1970),  and  Box 
and  Draper  (1971)  pointed  out,  however,  that  the  number  of  design  points  needed  to  adequately 
approximate  a  D-optimal  design  measure  might  be  too  large  as  compared  to  the  number  of  pa¬ 
rameters  in  the  model.  It  is  for  this  reason  that  discrete  D-optimal  designs,  that  were  mentioned 
earlier,  became  attractive.  We  refer  to  an  ;V- point  D-optimal  design  as  D„  -optimal  (see  St.  John 
and  Draper  1975).  It  is  important  here  to  point  out  that  the  Equivalence  Theorem  mentioned 
earlier  does  not  apply  to  Dy  -optimal  designs. 

Several  algorithms  are  available  to  construct  D„  -optimal  designs.  The  DETMAX  algorithm 
by  Mitchell  (1974)  and  the  one  by  Fedorov  (1972,  pp.  160-165)  generate  a  sequence  of  JV-point 
designs  with  nondecreasing  values  of  1  AT  A"  I  for  a  given  N.  Mitchell  and  Bayne  (1978)  constructed 
D-optimal  fractions  of  three-level  factorial  designs  for  second  order  models.  Welch  (1982)  presented 
a  catalog  of  DN  -optimal  designs  which  can  also  perform  well  by  other  design  criteria  for  a  small  loss 
in  the  D-optimality  efficiency.  A  review  of  algorithms  for  constructing  DH  -optimal  designs  was 
given  by  Cook  and  Nachtsheim  (1980). 

The  D*optimality  criterion  can  also  be  used  in  situations  where  a  subset  of  q  ( <  p)  parame¬ 
ters  in  a  model  are  of  interest,  whereas  the  remaining  p  -  q  parameters  are  regarded  as  nuisance 
parameters.  A  design  measure,  s*,  is  Z?,-optimal  if  it  minimizes  the  determinant  of  the  submatrix 
of  the  inverse  of  the  information  matrix  M(s)  which  corresponds  to  the  vector  of  the  q  parameters 
of  interest.  An  equivalence  theorem  similar  to  the  Kiefer  and  Wolfowitz  Equivalence  Theorem 
holds  for  Z>,-optimality  (see  Kiefer  1961). 
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A  lesser  known  optimality  criterion  is  the  c*optimality  criterion.  A  design  measure  j*  is  c* 
optimal  if  it  minimizes  the  variance  of  the  least-squares  estimator  of  a  linear  function  of  the  pa¬ 
rameter  vector  0.  Thus,  if  M{s)  is  the  information  matrix  and  £  is  a  given  vector  of  known  elements, 
then  s*  is  c-optimal  if  it  minimizes  c'Xt ~  '{s)c  over  the  class  of  all  design  measures.  Murty  and 
Studden  (1972)  used  this  criterion  to  obtain  optimal  designs  for  estimating  the  slope  of  a  response 
surface  at  a  given  point  of  the  experimental  region. 

Design  Efficiency 

Any  of  the  previously  discussed  optimal  designs  can  be  used  to  evaluate  and  compare  designs 
on  the  basis  of  efficiency.  For  example,  for  the  D-  and  G-optimality  criteria  the  corresponding  ef¬ 
ficiencies,  as  defined  by  Atwood  (1969),  are,  respectively, 

D-efficiency  -  pAf(r)|/sup  lAf(r)ljl/' 

(/-efficiency  =  pi  sup  v{£,  s), 
i«X 

where  s  is  a  given  design  measure,  p  is  the  number  of  parameters  in  the  model,  H  is  the  class  of  all 
design  measures,  and  v(s,  r)  is  the  standardized  prediction  variance  described  in  (2.8).  Note  that 
by  (2.12)  the  G-efficiency  (as  well  as  the  D-efficiency)  of  s  is  less  than  or  equal  to  unity  and  greater 
than  or  equal  to  zero.  Lucas  (1974,  1976,  1977)  used  the  D-  and  G-efficiencies  to  evaluate  the 
performance  of  some  traditional  response  surface  designs  for  second  order  models. 

Additional  Remarks  Concerning  Optimal  Designs 

The  various  optimality  criteria  described  earlier  are  usually  referred  to  as  alphabetic  design 
optimality.  This  term  was  originally  coined  by  Box  (1982)  who  expressed  certain  reservations  about 
'the  usefulness  of  this  approach  as  far  as  response  surface  designs  are  concerned*.  Researchers  in 
this  theory  have  produced  outstanding  results  based  on  intricate  mathematical  formulations,  but 
have  neglected  to  address  important  practical  difficulties  that  are  present  in  RSM  applications.  As 
a  result,  there  has  been  some  confusion  on  the  part  of  RSM  practitioners  who  attempt  to  use  the 
alphabetic  optimality  results. 


In  the  design  optimality  approach,  the  design  of  experiment  is  set  within  a  rigid  framework 
governed  by  a  set  of  assumptions,  the  most  prominent  of  which  eharactcrizes  the  fitted  model  as 
the  true  model  over  the  experimental  region.  It  is  rarely  the  case  when  an  experimenter  would  be 
willing  to  make  such  an  assumption.  Furthermore,  the  experimental  conditions  can  vary  in  time 
and  may  not  conform  to  assumptions  made  a  priori  at  the  outset  of  the  experiment.  This  has 
prompted  Box  et  al.  (1978,  p.  472)  to  state  that  'in  recent  years  the  study  of  optimal  design  has 
become  separated  from  real  experimentation  with  the  predictable  consequence  that  its  limitations 
have  not  been  stressed,  or,  often,  even  realized.' 

In  the  traditional  RSM  approach,  bias  suspected  of  being  present  in  the  fitted  model  plays  a 
significant  role.  This  is  the  basis  for  the  Box- Draper  criterion  discussed  in  Section  2.1.  In  fact  the 
concern  about  model  inadequacy  in  this  approach  has  prompted  other  researchers  to  consider  de* 
signs  solely  on  the  basis  of  bias  (see  Karson  et  al.  1969).  The  bias  criterion,  however,  provides 
protection  against  a  specific  form  of  model  inadequacy  according  to  what  we  postulate  as  the  true 
model.  Kiefer  (1975,  pp.  284-286)  criticized  certain  aspects  of  the  preoccupation  with  bias,  pointing 
out  examples  in  which  the  variance  criterion  is  compromised  for  the  sake  of  the  bias  criterion.  Both 
Kiefer  and  Box  agree  that  design  selection  should  be  guided  by  more  than  one  criterion,  'no  single 
simple  prescription  can  be  expected  to  yield  satisfactory  designs  in  all,  or  even  most,  applications' 
(Kiefer  1975,  p.  286)  (see  also  Box  1982,  Section  7). 

The  work  of  Fedorov  represents  a  'middle  ground'  between  Box's  and  Kiefer's  approaches. 
It  provides  the  mathematical  apparatus  for  the  derivation  of  the  theoretical  results,  as  well  as  the 
necessary  algorithms  needed  to  apply  them  so  that  they  are  accessible  to  comparatively  "broad  cir¬ 
cles  of  researchers  and  technologists'  as  Fedorov  (1972)  stated  in  the  preface  of  his  book  on  optimal 
experiments. 

In  addition  to  Fedorov's  book,  the  book  by  Silvey  (1980)  and  the  review  articles  by  Bandemer 
(1980)  and  Atkinson  (1982)  provide  recent  surveys  of  optimal  experimental  designs. 


2.6  Design  of  Mixture  Experiments 


Throughout  this  review  we  restate  our  frustration  that  much  of  the  technology  advanced  in 
RSM  seems  to  not  have  been  communicated  to  subject  matter  users  and  thus  they  are  'far  behind' 
recent  advances.  However,  the  area  of  mixture  experiments  is  an  exception.  From  the  early  pio¬ 
neering  work  in  the  late  50's,  techniques  in  the  use  of  mixture  designs  became  important  to  users 
of  RSM  where  the  response  to  a  mixture  of  ingredients  depends  on  the  relative  proportions  of  the 
ingredients.  These  experiments  are  characterized  as  follows:  If  x„  i  »  1,  2,  ... ,  k  represents  the 
proportion  of  the  rth  component  in  a  k  factor  combination  then  for  each  experimental  run 
0  £  xt  £  1,  i  -  1,  2,  ... ,  k  and  I  x,  -  1.  The  development  of  designs  and  the  analysis  of  data 
collected  under  these  constraints  has  been  the  subject  of  much  research.  Indeed,  the  particular  area 
has  been  the  subject  of  two  reviews  by  Cornell  (1973,  1979),  mentioned  prominently  in  two  addi¬ 
tional  reviews  by  Mead  and  Pike  (1975),  Steinberg  and  Hunter  (1984),  and  the  subject  of  a  text 
Cornell  (1981). 

/ 
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Some  highlights  in  the  development  of  this  area  follow.  Scheff&eaov  (1958)  is  credited  by 
most  with  having  initially  considered  the  mixture  problem  although  Claringbold  (1955)  recognized 
that  designs  for  studying  such  relationships  are  on  a  simplex.  Scbeff&eac.  proposed  the  simplex 
lattice  design  as  well  as  a  family  of  polynomial  models  to  use  in  approximating  the  underlying  re¬ 
lationship  between  the  response  and  mixture  components.  A  deficiency  of  the  simplex  lattice  design 

is  that  its  use  tends  to  result  in  mixtures  involving  only  two  components  regardless  of  the  total 

e' 

number  of  components  under  consideration.  To  remedy  this  Scheff&w.  (1963)  proposed  the 
simplex  centroid  design  which  consists  of  2*  -  1  design  points  which  are  the  centroid  of  the  simplex 
and  the  centroid  of  all  lower  dimensional  simplices  contained  within  the  original  (k  -  1)  dimen¬ 
sional  simplex. 

Thompson  and  Myers  (1968)  developed  designs  that  allow  the  use  of  the  minimum  average 
bias  and  variance  criteria  of  Box  and  Draper.  They  make  use  of  a  transformation  to  {k  -  1)  di¬ 
mensional  variables  so  that  the  standard  designs  can  be  used.  A  transformation  is  given  to  the  user 
which  produces  the  design  in  the  natural,  mixture  variables. 


It  is  not  unusual  for  a  mixture  experiment  to  be  such  that  there  are  bounds  on  some  or  all 
of  the  components  of  the  combination.  McLean  and  Anderson  (1966)  proposed  the  extreme 
vertices  design  for  such  situations.  Snee  (1981)  uses  the  gasoline  blending  problem  to  illustrate 
some  practical  problems  that  occur  when  considering  mixture  experiments.  Attention  is  focused 
on  the  design  of  the  blending  studies  and  appropriate  models  for  analysis  as  well  as  the  use  to  which 
the  fitted  model  is  put.  In  this  paper  currently  available  algorithms  are  used  in  the  construction  of 
linear  and  quadratic  model  designs  where  there  are  single  and  multiple  component  constraints  on 
the  region  of  feasible  blends.  Piepcl  (1983a)  discusses  several  guidelines  for  developing  constraints 
on  the  levels  of  components  in  mixtures  and  presents  techniques  for  checking  the  consistency  of  the 
constraints  that  are  developed.  Piepel  (1983b)  notes  that  the  centroid  of  a  constraint  region  in 
mixture  experiments  has  been  defined  as  the  average  of  all  extreme  vertices  of  the  region.  Utilizing 
the  definition  of  the  centroid  as  the  center  of  mass  or  volume  of  the  region  Piepel  develops  an  al¬ 
gorithm  to  calculate  centroids  of  various  dimensional  faces  of  the  constraint  region  as  well  as  the 
centroid  of  the  ( k  ~  1)  dimensional  simplex.  Whenever  constraints,  in  the  form  of  upper  and  lower 
bounds,  are  placed  on  the  components  of  a  mixture  the  size  and  possibly  the  slope  of  the  exper¬ 
imental  region  is  altered.  Crosier  (1984)  proposes  a  transforr  .ation  to  pseudo  components  to  re¬ 
duce  the  ill-conditioning  created  by  the  constraints.  Conditions  are  given  which  can  be  used  to 
determine  when  the  use  of  the  transformation  is  preferable  to  that  of  other  transformations. 

Hare  (1979)  considers  the  development  of  designs  for  the  situation  where  some  of  the  design 
variables,  called  process  variables,  are  not  subject  to  the  simplex  constraint,  1.  Of  concern 

here  is  the  choice  of  settings  for  the  mixture  variables  as  well  as  the  process  variables.  The  problem 
has  also  been  considered  by  Vuchkov  et  al.  (1981)  who  proposed  the  use  of  a  sequential  procedure 
to  produce  designs  that  were  nearly  D-optimal.  The  traditional  ScheffJbafc-  model  in  mixture  ex¬ 
periments  that  contain  process  variables  contain  cross-product  terms  between  the  mixture  and 
process  variables  as  well  as  indicate  that  the  process  by  mixture  cross-product  terms  estimate  the 
effects  of  the  process  variables  on  the  blending  properties  of  the  mixture  variables  only  and  do  not 
give  an  overall  measure  of  the  main  effects  and  interactions  of  the  process  variables  alone.  They 
develop  a  reparameterized  model  that  permits  the  experimenter  to  separate  the  effects  of  the  process 
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variables  from  the  mixture  variables.  From  such  a  model  it  is  demonstrated  that  reduced  model 
forms  can  be  obtained  through  the  use  of  a  variable  selection  procedure.- 

Aitchison  and  Bacon-Slone  (1984)  show  that  the  constant  sum  constraint,  i.e.  I  x,  »  1  can 

<•1 

be  removed  by  a  log  ratio  transformation,  z,  m  log— r-  /  -  1,  2,  ... ,  k  -  1.  When  this  is  done 
the  modeled  expected  response  can  be  expressed  as  a  polynomial  in  the  z,.  Since  the  z,  can  be  varied 
independently  the  polynomial  can  include  all  terms  of  appropriate  degree.  For  this  paper  attention 
is  restricted  to  linear  and  complete  quadratic  models. 

Rotalability  in  mixture  experiments  has  received  very  little  attention  by  comparison  to  other 
more  general  response  surface  investigations.  In  an  attempt  to  obtain  constant  prediction  variance 
within  the  simplex  factor  space,  Thompson  and  Myers  (1968)  developed  rotatable  response  surface 
designs  on  the  basis  of  their  transformation  of  the  k  mixture  variables  to  k  —  1  mathematically  in* 
dependent  variables.  These  designs  were  used  to  fit  polynomial  models  over  some  ellipsoidal  region 
inside  the  simplex  factor  space.  Cornell  and  Khuri  (1979)  introduced  a  different  configuration  of 
rotalability,  one  that  is  suited  for  ternary  mixture  problems.  In  this  new  configuration,  the  pre¬ 
dicted  response  has  constant  variance  on  each  triangle  of  a  set  of  concentric  triangles  within  the 
simplex  factor  space.  This  is  achieved  by  using  a  nonlinear  transformation  which  maps  concentric 
circles  in  a  new  coordinate  system  onto  concentric  equilateral  triangles  in  the  mixture  space. 
Standard  routable  response  surface  designs  in  the  new  coordinate  system  can  then  be  mapped  onto 
designs  in  the  mixture  space  resulting  in  the  prediction  variance  being  constant  on  concentric  tri¬ 
angles. 

2.7  Designs  for  the  Slope 

In  many  applications  of  RSM,  good  estimation  oi  the  derivatives  of  the  response  function 
may  be  as  important  or  perhaps  more  important  than  estimation  of  mean  response.  Certainly,  the 
computation  of  a  stationary  point  in  a  second  order  analysis,  or  the  use  of  gradient  techniques,  e.g., 
steepest  ascent  or  ridge  analysis  depend  heavily  on  the  partial  derivatives  of  the  estimated  response 
function  with  respect  to  the  design  variables.  Since  designs  that  attain  certain  properties  in  y  (esti¬ 
mated  response)  do  not  enjoy  the  same  properties  for  the  estimated  derivatives  (slopes),  it  is  im- 
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portant  for  the  user  to  consider  experimental  designs  that  are  constructed  with  the  derivatives  in 
mind. 

Atkinson  (1970)  considered  designs  for  estimation  of  the  slope  at  a  fixed  point  with  the  re* 
sponse  function  being  of  order  one.  The  criterion  used  is  the  expected  mean  squared  error  for  a 
directional  derivative,  averaged  over  all  possible  directions.  Ott  and  Mendenhall  (1972)  dealt  with 
the  special  case  of  a  single  design  variable  and  a  second  order  model.  No  model  misspccification 
was  considered  and  the  primary  focus  was  on  the  properties  of  the  variance  of  the  estimated  slope. 
They  considered  optimum  spacing  of  the  levels  of  the  design  variables  that  gave  desirable  properties 
of  the  variance  {unction.  Murty  and  Studden  (1972)  considered  polynomial  regression  models  with 
the  criterion  being  the  variance  of  an  estimated  slope  at  a  fixed  point  and  averaged  over  an  interval. 

Hader  and  Park  (1978)  extended  the  notion  of  rotatability  to  cover  the  slope  for  the  case  of 
second  order  models.  They  developed  design  parameters  for  central  composite  designs  that  result 
in  'slope  rotatability,'  i.e.,  the  variance  of  the  estimated  derivatives  is  constant  for  all  points 
equidistant  from  the  design  center.  Tables  are  produced  of  designs  for  values  of  k  -  2,  3,  ....  8. 
Some  of  the  designs  featured  have  replicated  axial  points  rather  than  heavy  replication  in  the  design 
center. 

Myers  and  Lahoda  ( 1975)  extend  the  Box  and  Draper  integrated  mean  squared  error  criterion 
under  model  misspccification  to  cover  sets  of  parametric  functions  with  the  slopes  being  primary 
applications.  For  the  second  order  case  they  develop  designs  that  are  minimum  slope  variance, 
conditional  on  achieving  minimum  integrated  slope  bias  in  case  of  both  spherical  and  cuboidal  re* 
gions  of  interest.  Their  designs  are  central  composite  designs  with  dual  axial  levels  and  are  tabulated 
for  k  ■  2 , 3 . 7. 

Mukerjee  and  Heda  (1985)  develop  designs  associated  with  minimum  variance  of  the  esti* 
mated  slope  maximized  over  all  points  in  the  factor  space  for  second  and  third  polynomial  models 
over  a  spherical  region. 
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2.8  Design  for  Minimum  Dias  Estimation 

K arson  et  al.  (1969)  introduced  the  notion  of  minimum  bias  estimation  in  response  surface 
analysis.  The  motivation  was  derived  from  the  work  of  Bo*  and  Draper  in  which  it  was  made  clear 
that  consideration  of  bias  in  RSM  is  at  least  as  important  as  consideration  of  variance.  As  a  result, 
an  alternative  to  least  so.uares  was  introduced  for  minimizing  average  bias  at  the  estimation  stage. 

Given  that  minimum  bias  estimators  are  used,  Karson  (1970)  introduced  a  design  criterion 
to  provide  protection  against  certain  higher  order  model  terms.  This  criterion  results  in  conditions 
on  the  design  moments  being  met.  With  a  single  independent  variable,  conditions  are  given  for 
fitting  up  to  a  cubic  polynomial  and  protecting  against  a  possible  fiAh  degree  polynomial.  The 
author  also  considers  the  two  factor  case  for  fitting  a  first  order  polynomial  and  protecting  against 
the  presence  of  second  and  third  order  terms.  A  general  p  factor  result  for  the  linear  to  quadratic 
case  is  also  given.  Karson  and  Spruill  (197S)  examine  four  alternative  design  criteria  given  that 
minimum  bias  estimation  is  used  in  the  case  of  a  single  independent  variable.  The  alternative  design 
criteria  are  minimaxes  or  minimum  averages  of  the  squared  differences  between  the  expected  value 
of  the  stationary  point  found  by  fitting  the  approximating  polynomial  and  either  a  prior  guess  at 
the  location  of  the  stationary  point  or  the  location  of  the  stationary  point  as  determined  from  the 
true  model.  Evans  and  Manson  (1978)  utilize  the  design  flexibility  associated  with  the  use  of  min¬ 
imum  bias  estimation  to  construct  D-optimal,  V-optimal,  and  A-optimal  experimental  designs  in 
two  dimensions  via  a  simplex  search  procedure.  The  designs  are  obtained  for  both  square  and  cir¬ 
cular  regions  of  interest. 

In  a  related  work,  Kupper  and  Meydrech  (1973)  consider  the  use  of  an  estimator  of  the  form 
ATfi  for  the  parameter  of  the  fitted  polynomial.  They  give  conditions  under  which  it  is  possible  to 
determine  K's  providing  smaller  J  than  when  K  ■  /  for  any  choice  of  experimental  design.  The 
case  where  the  fitted  polynomial  is  of  degree  one  and  the  true  polynomial  is  of  degree  two  is  con¬ 
sidered.  These  same  authors  (1974)  compare  their  approach  to  that  of  Box  and  Draper  (1959)  as 
well  as  that  due  to  Karson  et  al.  (1969). 


2.9  RSM  Design  Robustness 


As  we  indicated  earlier  in  this  paper,  the  first  notion  of  RSM  design  robustness  appeared  in 
the  work  by  Box  and  Draper  (1959,  1963)  dealing  in  protection  against  model  undcrspcciTication. 
It  formalized,  in  an  RSM  setting,  the  notion  of  proper  placing  of  the  design  points  in  from  the 
boundary  of  the  region  of  interest  when  underspecification  is  a  concern.  While  specifics  in  the  Box 
and  Draper  work  dealt  with  spherical  regions,  Draper  and  Lawrence  (1965)  applied  the  Box  and 
Draper  approach  to  generate  designs  robust  to  model  inadequacies  in  the  case  of  cuboidal  regions. 

Aside  from  model  misspecification,  RSM  design  robustness  includes  the  following  categories 

i)  Robustness  to  outliers  in  the  data 

ii)  Robustness  to  errors  in  the  design  levels 

iii)  Robustness  to  nonnormality  in  the  model  errors 

iv)  Designs  for  extrapolation  under  conditions  of  model  misspecifications. 

We  feel  committed  to  extend  the  notion  of  RSM  to  include  F-tests  that  are  developed 
through  regression  analysis.  Papers  by  Jensen  et  al.  (1975)  and  Vuchkov  and  Solakov  (1980) 
demonstrate  that  experimental  design  has  a  substantial  impact  on  the  robustness  to  normality  of 
the  F*test  in  a  regression  or  RSM  model.  In  the  former  paper,  the  class  of  orthogonal  designs  was 
found  to  be  robust  in  the  case  of  a  first  order  response  surface  model.  In  the  latter  case,  a  similar 
result  was  demonstrated.  The  approach  was  taken  of  considering  the  property  of  'quadratic  bal¬ 
ance'  in  a  design.  This  property  is  one  making  the  hat  diagonals  equal,  where  the  hat  diagonals 
are  the  diagonal  elements  ol 


H  -  XWX)"'* 

Khuri  arid  Myers  (1981)  implement  the  robustness  criterion  of  Box  and  Watson  (1962)  for  first 
order  RSM  models.  Again,  the  robustness  here  is  a  resistance  of  tests  of  significance  to  failure  of 
the  distribution  of  errors  to  be  normal.  A  method  is  described  for  construction  of  a  robust  design 
assuming  a  fixed  number  of  experimental  runs.  Here,  of  course  A’  is  an  /V  *  p  model  matrix  given 
in  the  general  linear  model 
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One  very  intriguing,  albeit  predictable,  result  from  the  Vuchkov  and  Solakov  work  is  that  the  most 
desirable  prescription  in  terms  of  this  type  of  robustness  is  a  design  with  uniformly  replicated  trials 
and  equal  hat  diagonals.  Example  designs  are  the  2*  factorial  arrays  and  appropriate  fractions. 

From  a  point  of  view  of  motivation,  the  paper  by  Box  and  Draper  (1975)  would  appear  to 

have  been  very  timely.  Since  the  early  70's  considerable  attention  has  been  paid  to  the  detection 

of  outiicrs  and  the  accommodation  of  outliers  in  regression  analysis.  This  work  by  Box  and  Draper 

linked  the  awareness  of  outliers  to  notions  of  RSM  experimental  designs.  Design  properties  were 

sought  that  resulted  in  'minimal  impact'  of  outliers  in  the  resulting  data.  The  result  was  the  choice 
n 

of  designs  for  which  £  h$  is  minimized  where,  again,  the  are  diagonals  of  the  hat  matrix.  Again, 
in  the  first  order  case,  the  orthogonal  designs  are  optimal. 

Herzberg  and  Andrews  (1976)  deal  with  optimal  designs  under  nonoptimal  conditions  such 
as  missing  observations  and  presence  of  outliers  while  dealing  mostly  in  designs  that  arc  constrained 
to  be  economical.  They  introduce  a  measure  of  robustness  which  is  called  a  'probability  of 
breakdown',  given  by 

Pr{\XDX\  -  0} 

where  D  is  a  diagonal  matrix  with  iih  main  diagonal  element 


0  with  probability  a  (a) 

1  with  probability  1  —  a(&) 


and  aCs)  is  the  probability  of  losing  a  data  point.  When  I X DX\  -  0,  the  coefficients  cannot  all 
be  estimated.  In  addition,  they  consider  £(  I  XDiX\ •")  as  an  alternative  measure  of  robustness.  In 
Andrews  and  Herzberg  (1979)  they  introduce  the  notion  of  efficiency  by  working  with  the  quantity 

1  E\XD2X\'Ip 
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where  M  is  the  Fisher  information  matrix  for  the  appropriate  D  optimal  design.  Examples  arc 
given  to  indicate  how  notions  like  this  can  be  used  to  select  robust  designs.  The  numerator  of  the 
above  expression  is  compared  for  central  composite  designs  with  varying  number  of  center  runs. 
The  probability  value  a(£)  is  varied  in  the  study  and  the  results  indicate  clearly  that  criteria  like 
these  can  be  used  to  provide  a  choice  between  designs  when  the  user  fears  a  loss  of  information  in 
the  form  of  missing  data. 

In  a  paper  that  considers  protection  against  model  misspecuication  and  outliers  simultane¬ 
ously,  Draper  and  Herzbcrg  (1979)  deal  with  n  integrated  mean  squared  error  criterion  much  like 
that  of  Box  and  Draper  in  1959  and  1963.  The  purpose,  though,  seemed  to  be  to  determine  if  the 
presence  of  outliers  produced  any  substantial  change  in  the  recommendations  made  earlier,  when 
one  only  considers  model  misspecification.  Special  cases  of  a  single  outlier  and  two  or  more  outliers 
are  considered.  As  expected,  in  the  case  of  a  first  order  fitted  model  with  protection  against  the 
existence  of  second  order  terms,  an  orthogonal  design  is  recommended.  In  addition,  the  recom¬ 
mendation  for  the  second  moment  is  to  make  it  slightly  larger  than  the  'minimum  bias'  value  un¬ 
covered  in  1959.  Specifically,  a  value  roughly  10%  larger  is  recommended.  Simply  put,  this  implies 
that  if  outliers  arc  present,  the  variance  portion  of  the  integrated  prediction  mean  square  error  be¬ 
comes  somewhat  more  important.  Extensions  to  the  case  of  the  second  order  fitted  model  were 
made  with  similar  recommendations. 

A  related  notion  of  robustness  deals  with  errors  in  the  factor  levels.  Vuchkov  and  Boyadjieva 
(1983)  consider  this  problem  and  attempt  to  determine  design  families  that  are  robust.  The  inter¬ 
ested  reader  should  first  read  Box  (1963)  in  which  the  effect  of  errors  in  factor  levels  is  considered 
in  both  fust  and  second  order  models.  While  no  specific  design  criteria  are  considered  in  this  work 
by  Box,  it  is  an  excellent  account  of  what  is  the  extent  of  the  damage  when  errors  in  design  variables 
occur. 


In  the  paper  by  Vuchkov  and  Boyadjieva,  some  of  the  same  concepts  covered  in  Box  are  re¬ 
stated.  They  define  an  estimator  that  is  BLUE  under  the  conditions  of  errors  in  factor  levels,  and 
which  depends  upon  the  moments  of  the  factor  levels.  Assuming  this  particular  form  of  estimation 
(which  reduces  to  least  squares  if  the  model  errors  are  unbiased),  they  define  'criteria  for  robustness' 
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and  make  some  comparisons  among  the  Box*Bchnken  designs,  rotatable  central  composite  designs, 
Hartley’s  small  composite  design,  Lucas’  optimal  composite  designs,  and  others. 

Finally,  we  deal  in  a  type  of  design  property  that  some  may  not  generally  classify  as 
robustness,  namely  designs  that  arc  resistant  to  errors  in  extrapolation.  This  is  particularly  important 
in  RSM  work  since  a  response  surface  is  often  used  by  necessity  for  extrapolation  purposes.  Draper 
and  Hcrzberg  (1973)  discuss  designs  that  are  robust  in  this  sense  under  the  setting  that  the  fitted 
model  is  first  order  and  one’s  ability  to  extrapolate  is  influenced  by  the  existence  of  second  order 
terms  in  the  true  structure.  Their  work  is  confined  to  consideration  of  extrapolation  in  a  specified 
direction  with  a  spherical  region  of  interest. 

Draper  and  Hcrzberg  (1979)  considered  both  the  first  and  second  order  case  in  which  interest 
is  in  extrapolation  outside  a  sphere.  The  region  of  extrapolation  is  a  spherical  shell.  Again,  variance 
and  bias  due  to  model  misspecification  are  considered  and  are  integrated  over  th*.  region.  In  the  first 
order  case,  they  reduce  the  problem  to  one  of  selecting  the  pure  second  mome'nt.  Ail  odd  moments 
through  order  three  are  assumed  zero.  They  are  able  to  find  the  second  design  moment  :hat  mini* 
mizes  what  they  call  a  'proportional  mean  square  error.'  Designs  with  maximum  spread,  i.t point* 
placed  on  the  perimeter  of  the  region  of  interest,  are  found  to  be  optimal.  In  the  second  order  case, 
their  work  was  restricted  to  second  order  rotatable  designs  with  all  odd  moments  less  than  or  equal 
to  five  being  zero.  Emphasis  is  put  on  the  central  composite  design.  The  conclusions  indicate  that 
the  rotatable  central  composite  designs  that  are  most  robust  to  extrapolation  require  design  points 
that  are  expanded  to  the  edge  of  the  region  of  interest. 

In  this  section  and  in  many  others  in  this  review,  the  influence  of  the  Box  and  Draper  pro¬ 
tection  against  model  inadequacy  becomes  evident.  What  we  have  seen  occur  in  much  of  the  design 
robustness  work  is  a  blending  of  the  problem  of  model  inadequacy  with  other  deviations  from 
ideality. 

Though  the  foregoing  might  be  viewed  as  representing  the  bulk  of  the  'robust  design'  work 
involved  strictly  in  the  use  of  RSM,  there  certainly  are  other  important  pieces  of  work  that  deal  in 
design  robustness  in  regression  models.  Atkinson  (1972)  discusses  designs  that  are  best  for  detection 
of  model  inadequacy.  Atwood  (1971,  1975)  deals  in  robust  procedures  for  estimating  polynomial 


regression  and  methods  for  estimating  a  response  surface  when  the  number  of  parameters  are  un¬ 
certain. 


One  important  note  of  caution  should  be  made  to  any  user  or  potential  user  of  RSM.  De¬ 
signs  are  chosen  because  they  are  optimal  or  near  optimal  in  some  sense,  or  because  they  are  robust 
or  near  robust  in  some  sense.  But  a  design  may  indeed  be  near  optima]  in  one  sense  (say  D* 
optimality)  but  not  near  optimal  with  respect  to  another  criterion.  Obviously,  a  desirable 
robustness  property  of  an  RSM  design  is  one  in  which  the  design  enjoys  a  near  optimal  or  at  least 
■’steady'  performance  in  several  criteria  that  is  important  to  the  user.  Tor  example,  Galil  and  Kiefer 
(1977)  and  Kiefer  (1975)  show  variation  in  performance  of  designs  as  one  moves  across  criteria. 

2.10  Designs  for  Special  Goals 

Designs  for  Nonlinear  Models 

A  nonlinear  model  is  a  model  of  the  form 

y  ■  fCs.S)  +  e  ,  (2.16) 

where  £  «  (jc„  Xj,  ... ,  x„Y  is  a  vector  of  design  variables,  2  ■  (9„  9„  ... ,  9,)'  is  a  vector  of  un¬ 
known  parameters,  e  is  a  random  error,  and  f  (js,  2)  is  a  known  function  nonlinear  with  respect  to 
at  least  one  element  of  2-  Nonlinear  models  have  not  received  a  great  deal  of  artention  in  RSM, 
even  though  they  have  many  applications  in  several  areas,  particularly,  in  biological  and  chemical 
sciences. 

The  main  design  criterion  for  nonlinear  models  is  the  D-optimality  criterion,  which  actually 
applies  to  a  linearized  form  of  the  model  in  (2.16)  (Box  and  Lucas  1959).  More  specifically,  if 
f  (i,  2)  is  approximately  linear  in  fi  in  a  neighborhood  of  fi  =  2o .  then  a  design  D  is  chosen  so  that 
the  determinant  I  F'(D,  &,)F{D,  2,)!  is  maximum,  where  F(D ,  &,)  is  an  N  *  p  matrix  whose 
(/,  y)th  element  is 


f,  =  fjfa:  8)  I 
'J  5*j  '2-fl, 


i  •  1.  2 . N- 

j  =  1.  2 . p. 


(2.17) 
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In  (2.17),  ~ •  denotes  the  partial  derivative  of  f  with  respect  to  0.,  &  is  the  value  of  a  in  the  /th 
c0/ 

experimental  run,  and  N  is  the  number  of  experimental  runs.  An  optimal  design  obtained  under 
this  criterion  is  more  appropriately  called  a  locally  D*optimal  design,  a  terminology  used  by 
ChemofT  (1953).  A  similar  design  criterion  is  available  if  the  interest  is  in  only  a  subset  of  the  pa* 
rameter  vector  fi  (see  Box  1971,  Hill  and  Hunter  1974). 

Unlike  D-optimal  designs  for  linear  models,  those  for  nonlinear  models  depend  upon  the 
unknown  values  of  2-  This  is  an  unappealing  characteristic  of  nonlinear  models  and  was  most  ap* 
propriately  depicted  by  Cochran  (1973):  "You  tell  me  the  value  of  fi  and  I  promise  to  design  the 
best  experiment  for  estimating  2.'  There  are  several  procedures  to  remedy,  or  at  least  alleviate,  the 
problem  of  design  dependency.  One  procedure  treats  2o  as  an  initial  guess  which  is  used  to  obtain 
a  point  design.  Observations  collected  at  the  points  of  this  design  are  then  used  to  provide  an 
estimate  of  2  •  Thereafter,  additional  design  points  are  augmented  sequentially  while  updating  the 
estimate  of  2  (see  Box  and  Hunter  1965).  On^  advantage  of  this  sequential  approach  is  that  it  helps 
reduce  the  dimensionality  of  the  design  optimization  problem.  If  the  initial  value,  2o.  of  2  is  a  good 
"guess',  then  under  certain  conditions,  Atkinson  and  Hunter  (1968)  showed  that  when  the  number 
of  design  points  N  is  a  multiple  of  p,  the  number  of  parameters,  the  optimal  JV-point  design  consists 
of  replications  of  the  optimal  p-point  design  obtained  on  the  basis  of  2o- 

An  alternative  strategy  for  dealing  with  design  dependency  is  to  adopt  a  Bayesian  approach. 
Zacks  (1977)  considered  maximizing  the  expected  value  of  the  determinant  I  F'(D,  2)  F(D,  2)1 
with  respect  to  some  prior  distribution  of  2-  Bayes  sequential  designs  can  also  be  obtained.  An 
optimal  design  for  a  new  stage  of  experimentation  is  determined  by  maximizing  the  aforementioned 
expected  value  with  respect  to  the  posterior  distribution  of  2  given  the  results  from  the  previous 
stages.  Another  strategy  that  is  applicable  when  the  number  of  design  variables  is  equal  to  one  is 
to  approximate  the  nonlinear  function  in  (2.16)  with  a  Lagrange  interpolating  polynomial  (Khuri 
1982).  The  interpolation  points  are  chosen  as  the  zeros  of  a  Chebyshev  polynomial  of  the  first  kind. 
Another  approximation  using  spline  functions  was  considered  by  Bumrungsup  (1984). 

The  dependency  of  the  optimal  nonlinear  design  on  the  parameter  vector  2  is  less  severe  if  the 
model  is  partially  nonlinear.  By  that  we  mean  that  some  elements  of  2  appear  linearly  in  the  model 
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while  others  appear  nonlinearly.  Hill  (1980)  showed  that  a  locally  D -optimal  design  for  a  partially 
nonlinear  model  depends  only  on  the  model's  nonlinear  parameters.  This  property,  however,  docs 
not  always  hold  when  the  design  is  for  estimating  a  subset  of  the  parameter  vector  (see  Khun  1584). 

The  use  of  the  D-optimality  criterion  for  choosing  a  nonlinear  design  may  not  always  be  ap¬ 
propriate  particularly  in  an  RSM  «etung.  One  must  remember  that  it  is  a  variance  criterion  which 
does  not  account  for  any  possible  bias  in  the  model.  Just  like  in  linear  models,  bias  is  an  important 
consideration  for  the  choice  of  design.  In  fact,  in  nonlinear  models  more  emphasis  should  be  given 
to  the  subject  of  bias  since  it  is  not  quite  clear  how  bias  should  be  defined.  Unfortunately,  little 
appears  to  be  known  about  this  subject.  A  brief  mention  of  bias  in  nonlinear  models  is  given  in  a 
paper  by  Atkinson  (1972,  Section  5). 

Designs  for  Fitting  Spline  Functions 

Let  (a,  6)  be  an  interval  that  is  partitioned  into  A  +  1  subintervals  by  the  points 
a  m  t,  <  t,  <  x,  <  —  <  ■  b.  A  spline  function  in  one  design  variable  x  is  a  continuous 

function  f(x)  defined  on  [a,  b\  such  that  f(jr)  restricted  to  [x,_„  t,],  i  »  1,  2,  ... ,  h  +  1  ,  is  a 
polynomial,  and  satisfies  certain  differentiability  conditions  at  the  points  t,,  t,,  ... ,  xA.  The 
polynomials  defined  on  these  subintervals  are  different,  and  the  partition  points, 
t„  t„  ....  and  x*  are  called  knots.  Spline  functions  are,  therefore,  segmented,  or  piecewise 
polynomials.  This  definition  of  spline  functions  can  be  extended  to  more  than  one  dimension. 

Spline  functions  are  useful  in  approximating  continuous  functions,  just  like  polynomials.  In 
certain  situations,  however,  they  are  more  suitable  than  polynomials.  This  is  particularly  true  when 
the  function  to  be  approximated  is  defined  over  a  wide  region  where  its  behavior  in  one  part  of  the 
region  is  unrelated  to  (or  undetermined  by)  its  behavior  in  another  part.  Polynomial  approximation 
on  the  other  hand  forces  some  kind  of  unwarranted  semblance  on  the  behavior  of  the  function. 
This  is  because  the  behavior  of  a  polynomial  over  an  entire  region  can  be  fully  determined  by  its 
behavior  over  only  a  part  of  the  region. 

Even  though  spline  functions  have  received  considerable  attention  in  the  mathematical  liter¬ 
ature,  their  introduction  into  the  statistical  literature  is  relatively  recent.  Wold  (1974)  gave  a  review 


of  the  use  of  spline  functions  in  least  squares  fitting  and  provided  several  useful  references  on  the 
subject.  Gallant  and  Puller  (1973)  discussed  the  estimation  of  the  location  of  the  knots  when  fitting 
a  spline  function  to  a  data  set.  Very  little  is  known  about  the  choice  of  design  for  fitting  a  spline 
function.  Studden  (1971)  discussed  D-optimal  and  L-optimal  designs  for  spline  functions.  Draper 
et  al.  (1977)  applied  the  Box* Draper  (1959)  design  criterion  for  protection  against  model  bias.  Park 
(1978)  adopted  the  D*optimality  criterion  when  fitting  a  spline  function  with  only  one  knot.  In  all 
of  these  papers,  the  knots  of  the  fitted  spline  function  were  considered  known.  More  recently, 
Bumrungsup  (1984)  used  spline  functions  to  approximate  the  mean  response  in  a  nonlinear  model. 
The  knots  were  chosen  so  that  the  error  of  approximation  did  not  exceed  some  specified  value 
chosen  a  priori.  The  approximating  spline  function  was  then  used  to  obtain  optimal  designs  for  the 
nonlinear  model.  These  designs  have  the  advantage  of  being  not  dependent  on  the  model's  pa¬ 
rameters. 

More  research  is  needed  in  this  area.  Optimal  designs  chosen  on  the  basis  of  spline  approxi¬ 
mation  depend  on  the  number  and  location  of  the  knots.  Even  if  the  number  of  knots  is  deter¬ 
mined  so  that  a  certain  degree  of  approximate  n  is  achieved,  the  positioning  of  the  knots  is  an 
important  consideration  which  should  be  addressed.  Extensions  to  more  than  one  dimension  are 
also  needed. 

Design  Augmentation 

Design  augmentation  is  an  increasingly  important  component  in  the  progression  of  RSM  re¬ 
search.  It  is  consistent  with  the  'sequential  spirit'  of  RSM.  In  a  large  percentage  of  RSM  appli¬ 
cations  more  than  one  stage  of  experimentation  and  analysis  is  required.  The  notion  of  moving 
from  one  experimental  region  to  another  and  planning  an  experiment  in  sequence  often  becomes 
necessary  in  real  RSM  situations.  Quite  often  additional  experimental  runs  are  taken  after  an  ex¬ 
periment  has  been  carried  out.  This  may  be  necessary  in  the  following  situations: 

i)  When  the  initial  experiment  is  poorly  planned  and,  consequently,  results  in  observa¬ 
tions  that  are  not  as  informative  about  the  fitted  model  as  is  desired. 


ii)  When  the  prescribed  design  settings  are  not  attained  by  the  actual  values  of  the  design 
variables  used  in  the  initial  experiment.  This  might  be  caused,  for  example,  by  tech* 
nic&l  limitation  or  by  equipment  malfunction. 

iii)  When  there  is  a  need  to  fit  a  model  that  is  more  complex  than  the  original  model, 
which  may  be  inadequate,  and  the  initial  design  does  not  accommodate  such  an  action. 

iv)  When  the  experimenter  decides  to  change,  or  modify,  the  region  of  interest  in  the  ori¬ 
ginal  experiment.  This  usually  occurs  in  sequential  experimentation  where  information 
obtained  in  one  stage  is  used  to  plan  the  next  stage.  Also,  there  might  be  a  need  to 
gain  more  information  about  the  response  in  certain  parts  of  the  region  that  are  not 
adequately  covered  by  the  initial  design. 

In  any  of  the  above  situations,  the  additional  runs  are  chosen  in  a  manner  that  remedies,  or 
repairs,  the  initial  design.  Situation  (i)  was  considered  by  several  authors.  Dykstra  (1966)  and 
Gaylor  and  Merrill  (1968)  augmented  a  nonorthogonal  design  for  fitting  a  first  order  model  with 
addi'ional  runs  in  order  to  reduce  multicollinearity  effects  among  the  design  variables.  Dykstra 
(1966)  suggested  several  other  criteria  for  augmenting  a  nonorthogonal  first  order  design.  Covey- 
Crump  and  Silvey  (1970)  and  Mayer  and  Hendrickson  (1973)  adopted  the  maximization  of  \X"X\ 
criterion  for  the  selection  of  the  additional  experimental  runs,  where  X  is  the  matrix  of  design  var¬ 
iables  associated  with  the  linear  model  for  both  the  initial  design  and  the  augmented  runs.  In 
Covey-Crump  and  Silvey  (1970)  the  augmented  design  points  were  selected  from  a  spherical  region, 
whereas  in  Mayer  and  Hendrickson  (1973),  the  augmented  design  points  were  constrained  by  cost 
considerations.  Dykstra  (1971)  described  a  method  for  augmenting  the  initial  design  with  additional 
points  taken  one  at  a  time.  Each  point  is  chosen  from  a  list  of  candidates,  the  selection  criterion 
being  the  choice  of  the  location  for  which  prediction  variance  is  maximized.  This  is  equivalent  to 
the  maximization  of  \X"X\  with  each  added  point.  Hebble  and  Mitchell  (1972)  followed  a  similar 
procedure,  except  that  the  entire  region  of  interest,  instead  of  a  set  of  candidate  points,  is  searched 
for  the  next  design  point.  Wynn  (1970)  showed  that  this  process  leads  in  the  limit  to  a  D-optimal 
design.  Hebble  and  Mitchell's  (1972)  procedure  can  also  be  useful  in  situations  (iii)  and  (iv).  Evans 
(1979)  used  the  simplex  search  technique  to  add  several  points  simultaneously  under  the  maximum 
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of  I  A"  A' I  criterion.  Suich  and  Derringer  (1977),  in  their  discussion  of  the  adequacy  of  the  regression 
model,  considered  data  augmentation  as  a  means  to  satisfy  their  criterion  of  adequacy.  They  did 
not,  however,  provide  a  systematic  procedure  for  the  choice  of  the  added  points. 

The  main  criterion  used  in  design  augmentation  has  been  the  maximization  of  lATAfl.  Re¬ 
cently,  Khuri  (1985a)  introduced  a  measure  of  rotatability  which  can  be  used  to  quantitatively  as¬ 
sess  departure  from  rotatability  of  a  nonrotatable  design.  Using  rotatability  as  a  criterion,  this 
measure  can  be  utilized  to.  repair  a  nonrotatable  design  by  design  augmentation.  Hebble  and 
Mitchell  (1972)  have  indirectly  restored  rotatability  through  their  maximization  of  lATAl  over  a 
spherical  region.  The  measure  of  rotatability  is  particularly  useful  in  situation  (ii)  where  a  design, 
originally  planned  to  be  rotatable,  becomes  nonrotatable  due  to  failure  to  attain  the  prescribed  de¬ 
sign  settings  for  a  rotatable  design.  It  can  also  be  used  in  situation  (iv)  where  a  rotatable  design  is 
intentionally  modified  by  design  augmentation  in  order  to  gain  more  information  in  certain  areas 
of  interest  (see  Littell  and  Mott  1974).  Such  a  modification  usually  results  in  loss  of  rotatability. 

Designs  to  Increase  the  Power  of  the  Lack  of  Fit  Test 

Checking  the  adequacy  of  a  fitted  model  is  an  important  consideration  in  RSM.  Draper  and 
Herzberg  (1971)  provided  some  insight  into  the  nature  of  hypotheses  that  can  be  tested  by  the  lack 
of  fit  test  when  replicate  observations  are  available  at  one  or  more  design  points.  Equally  important 
is  the  ability  to  detect  lack  of  fit  in  a  model  at  an  early  stage  of  the  experimental  process.  This 
necessitates  the  use  of  a  design  which  can  induce  a  certain  degree  of  sensitivity  in  the  lack  of  fit  test. 
To  accomplish  this  the  design  should  be  selected  so  as  to  maximize  the  power  of  the  associated  lack 
of  fit  test. 

Suppose  that  at  some  stage  of  the  experiment  the  following  model  is  Gtted: 

U(i)  ■  *,fii  .  (2.18) 

where  qOc)  denotes  the  mean  response  vector.  The  random  cnors  associated  with  the  observed 
response  values  are  assumed  to  be  independently  and  identically  nom  illy  distributed  with  zero 
mean  and  constant  variance.  The  true  model  is  assumed  to  have  the  form 
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The  noncentiality  parameter  associated  with  the  lack  of  fit  test  for  model  (2. 1 8)  is  proportional  to 
the  quadratic  form 

9  -  ,  (2-20) 


where 


In  order  to  maximize  the  power  of  the  lack  of  fit  test,  the  design  should,  therefore,  be  cboscn  to 
make  0  as  large  as  possible. 

Unfortunately,  designs  chosen  under  the  above  criterion  depend  on  &  ,  which  is  unknown. 
Several  approaches  were  considered  to  overcome  this  problem.  Atkinson  (1972)  proposed  the 
maximization  of  the  determinant  I L I .  This  is  the  same  as  the  D, -optimality  criterion  discussed  in 
Section  2.5  for  the  precise  estimation  of  a  subset  of  the  parameters  consisting  of  the  elements  of 
Designs  chosen  in  this  manner,  however,  may  provide  poor  estimates  of  fi,  if  the  fitted  model 
is  adequate.  Atkinson  (1972)  suggested  the  maximization  of  14 1  subject  to  maintaining  «  lower 
bound  on  the  efficiency  of  the  lack  of  fit  design  for  estimating  .  Atkinson  and  Fedorov's  (1975) 
T-optimality  criterion  is  essentially  the  same  as  the  one  we  have  mentioned.  They  suggested  three 
possible  solutions  to  deal  with  the  problem  of  dependency  on  ;  these  are  (i)  a  sequential  proce¬ 
dure  to  construct  a  T-optimal  design  whereby  observations  obtained  at  a  given  stage  are  used  to 
estimate  and  then  select  the  next  design  point,  (ii)  a  Bayesian  solution  which  requires  specifica¬ 
tion  of  a  prior  distribution  on  £3 ,  (iii)  a  maximin  procedure  in  which  the  minimum  value  of  G  over 
a  convex  region  9  in  the  £;-space  is  maximized  by  a  proper  choice  of  design.  A  special  case  of  the 
maximin  procedure  was  adopted  by  Jones  and  Mitchell  (1978)  using  a  particular  convex  region  <p 
determined  by  the  moment  matrices  of  the  experimental  region.  They  called  a  design  obtained 
under  the  maximin  crit  irion  a  A, -optimal  design.  They  also  considered  another  variant  of  the 
maximin  procedure,  namely,  the  maximization  of  the  aveiage  of  0  over  the  boundary  of  the  region 
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<p.  A  design  obtained  under  this  criterion  is  called  a  Aj  -optimal  design.  The  latter  criterion  is 
prefened  over  the  A, -optimality  criterion  which  cannot  be  used  in  situations  where  the  minimum 
value  of  0  over  q>  is  equal  to  zero. 

The  Aj-optimality  criterion  was  also  applied  by  Morris  and  Mitchell  (1983)  to  obtain  designs 
to  detect  the  presence  of  two-factor  interactions  among  several  two-level  factors  when  fitting  a  first 
order  model  in  these  factors.  An  extension  of  the  A,-optimality  criterion  as  well  as  the  A, -opti¬ 
mality  criterion  to  a  multircsponse  situation  was  recently  developed  by  Wijesinha  and  Khuri 
7 

(198»). 

The  traditional  lack  of  Gt  test  mentioned  earlier  requires  the  availability  of  replicate  observa¬ 
tions.  This  allows  the  partitioning  of  the  residual  sum  of  squares  from  the  analysis  of  the  Gtted 
model  into  a  sum  of  squares  due  to  lack  of  fit  and  a  sum  of  squares  due  to  pure  error.  When  such 
a  partitioning  is  not  feasible  (this  occurs  when  the  number  of  terms  in  the  Gtted  model  is  equal  to 
the  number  of  distinct  design  points),  the  method  of  check  points  (see  Schefle  1958,  Snee  1977)  can 
be  used  instead.  Check  points  ?re  points  in  the  experimental  region  where  data  are  collected  and 
compared  to  the  predicted  response  values  (at  these  points)  as  obtained  under  the  fitted  model. 
Shelton  et  al.  ( 1983)  introduced  a  method  for  selecting  the  check  points  so  as  to  maximize  the 
power  of  Scheffe's  (1958)  lack  of  fit  test  statistic,  which  makes  use  of  the  residuals  at  the  check 
points. 


III.  RESPONSE  SURFACE  ANALYSIS 

In  the  Hill  and  Hunter  review  of  1966,  four  steps  were  outlined  as  part  of  standard  procedure 
in  a  response  surface  analysis.  We  think  they  are  worth  repeating  here.  They  are  (i)  perform  a 
statistically  designed  experiment,  (ii)  estimate  the  coefficients  in  the  respoase  surface  equation,  (iii) 
check  on  the  adequacy  of  the  equation  (via  a  lack  of  Gt  test),  and  (iv)  study  the  response  surface  in 
the  region  of  interest.  It  is  (iv)  that  we  plan  to  consider  in  this  sect, on.  We  consider  subsections 
on  multiple  responses,  sampling  properties  of  optima,  and  models  that  are  alternatives  to  the 
standard  polynomial  models.  However,  underlying  all  of  these  topics  is  the  question  'What  do  we 
mean  by  a  response  surface  exploration?' 
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The  purpose  of  an  RSM  analysis  is  to  answer  certain  general  questions  regarding  the  nature 
of  the  response  function  in  the  region  of  the  experimental  design  Determining  whether  the  system 
contains  a  maximum  or  minimum,  or  is  a  saddle  system  may  be  very  important  to  the  exper¬ 
imenter.  On  the  other  hand,  the  existence  of  a  ridge  system  may  indicate  the  need  for  further  ex* 
perimentation  in  another  region.  In  some  situations,  response  surface  analysis  serves  as  more  than 
merely  an  empirical  device.  It  can  lead  to  what  Hill  and  Hunter  call  an  'elucidation  of  an  under* 
lying  mechanism',  and  thus  uncover  theoretical  information. 

A  distinction  should  be  made  between  analyses  whose  goals  are  to  explore  the  response  sux* 
face  as  indicated  above  and  that  in  which  the  major  goal  is  quality  estimation  of  optimum 
conditions.  The  method  of  steepest  ascent  remains  a  viable  technique  for  sequentially  moving  to* 
ward  the  optimum  response.  The  reader  is  referred  to  Davies  (1954)  and  Myen  (1976)  for  details 
and  examples.  Brooks  and  Mickey  (1961)  and  Myers  and  Khuri  (1979)  discuss  strategy  in  steepest 
ascent. 

It  is  interesting  to  note  that  in  the  Mead  and  Pike  review,  they  note  that  there  were  relatively 
few  applications  of  response  surface  optimization  with  the  use  of  procedures  other  than  steepest 
ascent.  At  this  point,  1 1  years  later,  we  might  echo  a  similar  observation  but  add  that  there  have 
been  applications  with  the  use  of  the  Nelder  Mead  (1965)  simplex  procedure.  What  still  surpris* 
ingly  appears  to  be  lacking  are  many  examples  of  formal  RSM  optimization  in  which  a  response 
model  other  than  first  order  is  used.  One  of  Mead  and  Pike's  conjectures  may  still  be  valid,  namely 
there  is  a  general  conservative  nature  in  many  research  workers.  Of  course,  the  theory  associated 
with  most  numerical  optimization  methods  is  developed  in  the  literature  under  an  assumed 
deterministic  model.  There  are  overriding  reasons  why  the  optimization  phase  of  RSM  is  not 
formally  applied  in  practical  situations.  There  is  a  general  unawareness  among  RSM  users  that 
optimization  techniques  are  available.  And  certainly  the  commercially  available  software  that  data 
analysts  use  for  other  phases  of  analysis  generally  does  not  contain  the  optimization  phase.  Nu¬ 
merical  analysis  has  advanced  to  the  extent  that  optimization  of  second  order  or  certain  nonlinear 
response  functions  with  or  without  constraints  presents  no  serious  difficulty.  And  yet  how  many 
software  packages  put  together  a  RSM  package  complete  with  canonical  analysis  and  optimization 
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procedures?  Indeed,  we  have  experienced  a  real  need  for  such  packages  in  industry.  We  attempt 
to  shed  more  light  on  this  subject  in  subsequent  sections. 

One  useful  procedure  in  a  second  order  response  surface  analysis  is  ridge  analysis,  introduced 
by  Hoerl  (1959)  and  formalized  by  Draper  (1963).  The  method  is  discussed  with  examples  in  Myers 
(1976).  Hoerl  (1985)  gives  a  historical,  almost  anecdotal  account  of  his  father's  development  of 
ridge  analysis.  Reasons  are  given  for  what  has  appeared  to  be  a  relative  obscurity  of  this  procedure, 
which  essentially  is  a  steepest  ascent  procedure  for  second  order  models.  It  serves  as  a  practical 
optimization  procedure  in  which  maximum  (or  minimum)  values  of  y  are  determined  on  radii 
varying  from  zero  to  the  perimeter  of  the  design  region.  The  resulting  locus  of  points  not  only 
produces  an  alternative  to  canonical  analysis  for  exploring  the  response  surface  but  it  also  can  offer 
a  recommended  set  of  operating  conditions  often  sought  by  the  analyst  when  the  canonical  analysis 
does  not  indicate  an  optimum  at  a  stationary  point  within  the  design  region.  With  the  entire  locus 
of  points  produced  by  the  maximum  ridge  and  a  knowledge  of  the  properties  of  the  experimental 
design  at  or  near  the  design  perimeter,  the  user  becomes  armed  with  information  that  allows  a 
reasonable  recommendation. 

It  has  been  our  observation  that,  in  practical  RSM  problems,  occasions  in  which  the  sta¬ 
tionary  point  produces  a  satisfactory  result  are  rare.  More  often  than  not,  ridge  conditions  or  saddle 
systems  are  experienced.  As  a  result,  ridge  analysis,  the  notion  of  finding  best  conditions  subject 
to  'anchoring'  the  results  inside  the  experimental  region,  or  some  other  type  of  constrained  opti¬ 
mization,  is  often  indicated.  The  advantage  in  ridge  analysis  as  an  exploratory  tool  has  intuitive 
appeal.  Most  users,  e.g.,  engineers,  biologists,  etc.,  understand  it  and  would  certainly  prefer  it  to 
drawing  conclusions  from  response  contour  plots  for  cases  where  k  i  3.  Hoerl  ( 1985)  laments  the 
lack  of  attention  to  ridge  analysis  in  the  statistics  literature.  He  walks  the  reader  through  a  nice 
practical  application.  Smith  (1976)  relates  situations  in  which  ridge  analysis  is  misleading.  Khuri 
and  Myers  (1979)  offer  modifications  to  ridge  analysis  in  cases  where  the  design  is  not  aear 
rotatable.  They  suggest  that  the  constraint  applied  should  be  on  the  prediction  variance  and  point 
out  that  if  the  design  is  not  rotatable,  the  modification  provides  more  desirable  results  for  the  user. 
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In  the  subjections  that  follow,  response  surface  analysis  in  the  ease  of  multiple  responses  is 
reviewed.  Doth  exploration  of  the  response  and  estimation  of  coefficients  are  considered.  This  is 
followed  by  the  vitally  important  topic  of  sampling  properties  of  estimated  optimum  conditions. 
Mead  and  Pike,  in  the  1975  review,  indicated  that  statistical  research  into  experimental  design  (say 
RSM  design)  has  outstripped  research  in  the  area  of  analysis.  To  a  great  extent,  this  remains  true. 

A  quick  count  of  the  volume  of  our  bibliography  in  design  as  opposed  to  papers  in  analysis  (in 
statistical  journals)  will  underscore  the  point. 

3. 1  Multiresponse  Analysis 

In  many  experimental  situations  a  number  of  responses  can  be  measured  at  a  time  for  each 
setting  of  a  group  of  design  variables.  Examples  of  multiresponse  situations  are  numerous.  In 
chemical  engineering  kinetics,  several  reactants  might  be  involved  in  a  chemical  reaction  which  give 
rise  to  a  number  of  responses  that  can  be  measured  simultaneously.  These  responses  are  usually 
represented  by  nonlinear  mechanistic  models,  or  by  linear  empirical  models  whenever  the  reaction 
mechanism  is  too  complex  such  as  in  an  industrial  reaction.  Ziegel  and  Gorman  (1980)  gave  an 
interesting  expose  of  the  use  of  multiresponse  data  for  studying  kinetic  models  in  the  petroleum 
industry.  The  area  of  foods  and  other  consumer  product  research  also  represents  a  challenge  be¬ 
cause  of  the  need  to  almost  always  study  many  responses.  Typical  studies  require  panel-type  re¬ 
sponses  or  'scores'  on  many  features  of  the  product.  The  review  by  Hill  and  Hunter  (1966)  listed 
several  examples  of  multiresponse  experiments  and  is  perhaps  the  first  paper  in  the  statistical  liter¬ 
ature  to  direct  attention  to  their  important  applications. 

Subsequent  papers  on  this  subject  have  stressed  the  utility  of  analyzing  multiresponse  data 
by  means  of  multivariate  techniques  which  take  into  account  interrelations  among  the  responses. 
A  univariate  analysis  in  which  responses  are  analyzed  individually  does  not  recognize  such  interre¬ 
lations  and  will  undoubtedly  fail  to  adequately  describe  the  true  mechanism  underlying  the  multi¬ 
response  system.  Unfortunately,  of  the  very  few  techniques  that  are  currently  available  in  the 
statistical  literature  for  the  analysis  of  multiresponse  data,  even  feweT  are  actually  used  by  practi¬ 
tioners.  Our  search  for  practical  RSM  applications  reveal  that  users  are  far  from  up-to-date  in  the 


40 


use  of  multircsponsc  techniques.  (See  Section  IV.)  It  is,  therefore,  imperative  that  existing  as  well 
as  future  multiresponsc  techniques  be  made  accessible  to  data  analysts.  It  is  also  equally  important 
to  demonstrate  the  utility  of  using  the  multiresponsc  approach  when  several  responses  are  analyzed. 

Multiresponse  Estimation 

The  general  multiresponse  model  is 

u  -  1,2 . N-, 

yut  *  G)  +  (3.i) 

i  -  1.2 . r 

where  a,  is  the  vector  of  settings  of  k  design  variables  at  the  uth  experimental  run,  B  is  a  vector  of 
unknown  parameters,  f,  is  a  function  of  known  form  for  the  rth  responses,  and  e„  is  a  random  error 
associated  with  the  /th  response  for  the  experimental  run  u.  It  is  assumed  that  the  e,/s  are  normally 

'  distributed  such  that  £(e„)  ■  0,  £(c,(  e,y)  •  0  for  all  /,  J,  u  #  v;  Var(E  J  »  au,  l  ■  1,2 . r\ 

£(e„  ct/)  *  a,  for  all  u,  i  *  j.  If  the  o,'s  are  known,  then  an  estimate  of  ft  can  be  obtained  by  using 
weighted  least  squares.  Most  often,  however,  the  o,'s  are  unknown.  In  this  case,  Box  and  Draper 
(1965)  using  the  Bayesian  approach  and  assuming  noninformative  prior  distributions  for  ft  and  the 
Oij  's  were  able  to  show  that  an  estimate  of  ft  can  be  obtained  by  minimizing  the  determinant 
1 5(fi)  1  with  respect  to  ft,  where 

S(&)  •  (Y  -  F)\Y  -  F)  ■  (3.2) 

In  (3.2),  Y  -  :£,]  is  the  data  matrix  with  %  being  the  vector  of  rth  response  values 

(/  ■  1,  2,  ... ,  r)  ,  and  F  is  an  N  *  r  matrix  whose  (u,  i)th  element  is  f,  (&,,  ft)  .  The  method  of 
estimation,  which  is  basically  equivalent  to  the  maximization  of  the  marginal  posterior  density  of 
ft  (under  the  assumed  noninformative  priors),  is  referred  to  as  the  Box-Draper  estimation  criterion. 
This  is  a  general  criterion  and  applies  to  models  whose  functional  forms  can  be  either  linear  or 
nonlinear  in  the  parameters. 

Box  et  al.  (1973)  cautioned  that  the  Box- Draper  estimation  criterion  can  lead  to  meaningless 
results  when  exact  linear  relationships  exist  among  the  responses.  They  called  such  relationships 
stoichiometric.  They  also  pointed  out  that  small  rounding  errors  in  the  responses  can  cause 
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15(0)1  to  be  different  from  zero  and  to  change  as  the  elements  of  0  arc  changed,  even  under 
stoichiometric  relationships.  To  resolve  this  confusing  situation,  they  introduced  an  eigenvalue 
analysis  which  checks  for  the  possibility  of  5(0)  having  a  zero  eigenvalue  aAcr  accounting  for 
round-ofT  errors.  They  wens  able  to  accomplish  this  by  devising  a  detection  procedure  which  in¬ 
volves  the  examination  of  the  eigenvalues  of  the  matrix  DD‘,  where  D  is  the  r  *  N  matrix 

d  -  r\i  -  it/ni .  (3.3) 

In  (3.3),  I  is  the  identity  matrix  of  order  N  *  N  and  1  is  the  vector  of  ones  of  order  N  *  1.  Box 
et  al.  ( 1973)  demonstrated  that  m  linearly  independent  relationships  exist  among  the  responses  if 
and  only  if  the  matrix  DD'  has  a  zero  eigenvalue  of  multiplicity  m.  The  orthonormal  eigenvectors 
of  DD'  which  correspond  to  a  zero  eigenvalue  of  DD '  are  used  to  identify  linear  relationships  among 
the  responses.  These  relationships  can  then  be  used  to  drop  m  responses  (if  m  linear  relationships 
exist  among  the  responses)  so  that  the  remaining  r  -  m  responses  are  not  linearly  related.  It  is  not 
always  easy  to  determine  which  m  responses  to  drop.  In  many  cases,  the  structure  of  the  problem 
may  well  dictate  a  natural  way  for  dropping  responses.  Having  selected  a  subset  oi  r  -  m  linearly 
independent  responses,  the  estimation  process  can  then  proceed  by  applying  the  Box- Draper  crite¬ 
rion  to  these  responses. 

Further  problems  associated  with  the  use  of  the  Box- Draper  estimation  criterion  were  con¬ 
sidered  by  McLean  et  al.  (1979).  Box  et  al.  (1970)  and  Stewart  and  Sorensen  (1981)  discussed  ap¬ 
plications  of  the  criterion  when  some  of  the  observations  are  missing.  More  recently,  Bates  and 
Watts  (1985)  proposed  a  new  computing  method  for  the  determination  of  the  parameter  estimates 
using  this  criterion.  The  method  is  based  on  a  generalized  Gauss* Newton  algorithm  for  the  min¬ 
imization  of  1 5(0)  I .  The  partial  derivatives  of  the  expected  responses  (with  respect  to  the  elements 
of  0),  which  are  required  in  this  algorithm,  are  generated  automatically  whenever  the  response 
models  are  determined  by  dynamic  models  defined  by  a  system  of  differential  eauations.  An  ap¬ 
proximate  confidence  region  on  the  parameter  vector  0  was  also  described  and  developed  by  Bates 
and  Watts  (1985)  using  a  quadratic  approximation  of  15(0)1 . 
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Design  of  Experiments  for  Multiresponse  Situations 

This  is  an  important  aspect  of  multircsponse  analysis,  yet  it  is  probably  the  least  developed. 
We  would  be  understating  by  indicating  that  the  development  of  multiresponse  designs  has  been 
lagging.  Yet,  as  we  reveal  in  subsequent  sections  on  practical  applications,  nearly  all  RSM  appli¬ 
cations  are  multiresponse  in  nature.  In  a  multiresponse  situation,  the  choice  of  design  should  be 
based  on  a  criterion  which  involves  all  the  responses. 

Draper  and  Hunter  (1966,  1967)  obtained  design  criteria  for  parameter  estimation  for  models 
of  the  form  described  in  (3.1).  These  criteria,  however,  require  knowledge  of  Z,  the  variance- 
covariance  matrix  for  the  r  responses.  Fedorov  (1972)  introduced  an  algorithm  for  the  construction 
of  a  D-optimal  design  for  a  linear  multiresponse  model  using  a  procedure  whereby  design  points 
are  chosen  sequentially.  Fedorov's  algorithm  also  requires  knowledge  of  I  .  Recently,  Wijcsinha 
and  Khuri  (1986)  introduced  a  modification  of  Fedorov's  algorithm  which  can  be  used  when  Z  is 
not  known. 

Another  criterion  for  the  choice  of  a  multiresponse  design  for  a  linear  multiresponse  model 
is  based  on  the  multivariate  lack  of  fit  test  developed  by  Khuri  (1985b).  The  design  is  chosen  so 
as  to  increase  the  power  of  this  test.  Details  concerning  the  construction  of  such  a  design  are  given 
in  Wijcsinha  and  Khuri  (1987). 

A  Test  for  Lack  of  Fit  of  a  Linear  Multiresponse  Model 

When  fitting  a  multiresponse  model,  provision  should  be  made  to  test  the  adequacy  of  the 
model  to  represent  the  behavior  of  the  multiresponse  system.  Since,  in  general,  the  responses  can 
be  correlated,  one  should  avoid  applying  the  usual  single-response  lack  of  fit  test  to  each  response 
model  individually.  Lack  of  fit  in  one  response  model  may  influence  the  fit  of  the  other  responses. 
Therefore,  multivariate  techniques  are  required  to  assess  the  overall  adequacy  of  the  response 
models  under  consideration.  Khuri  (1985b)  presented  such  a  test  for  linear  multiresponse  models. 


Multiresponse  Optimization 


When  several  responses  are  considered  simultaneously!  the  problem  of  determining  the  op* 
timum*  is  undefined  until  some  multivariate  optimization  criterion  has  been  chosen.  In  a  multire- 
sponse  situation,  no  unique  way  exists  for  ordering  values  of  a  multircsponse  function. 
Furthermore,  conditions  which  are  optimal  for  one  response  may  be  far  from  optimal  or  even 
physically  impractical  for  the  other  responses.  It  is  interesting  that  many  (perhaps  most)  users  re* 
quire  multiresponse  optimization  and  the  primary  medium  for  doing  the  analysis  is  simply  over¬ 
lapping  response  contour  maps.  Evidence  of  this  will  be  revealed  in  the  treatment  of  subject  matter 
applications  in  IV. 

Graphical  methods  were  used  early  in  the  history  of  RSM  development.  By  superimposing 
response  contours  and  visually  searching  for  a  common  region  where  the  responses  achieve  near 
optimal  values  (if  such  a  region  exists),  a  practitioner  might  be  able  to  somehow  arrive  at  a  location 
(or  locations)  of  a  'compromised'  optimum.  Lind  et  al.  (1960)  used  this- procedure  to  obtain  op¬ 
erating  conditions  for  maximizing  the  yield  and  minimizing  the  cost  of  a  certain  antibiotic.  Obvi¬ 
ously.  this  procedure  is  difficult,  if  not  impossible,  to  apply  when  the  number  of  design  variables 
exceeds  three.  Myers  and  Carter  (1973)  introduced  an  algorithm  for  determining  conditions  on  the 
design  variables  which  maximize  a  'primary  response'  function  subject  to  the  condition  that  a 
'secondary  response'  function  not  exceeds  a  certain  value.  The  development  of  this  algorithm  is 
somewhat  similar  to  that  of  ridge  analysis  used  with  a  single-response  model.  Biles  (1975)  extended 
this  concept  of  constrained  optimization  by  placing  bounds  on  the  values  of  several  'secondary  re¬ 
sponse'  functions.  Biles'  procedure  employs  a  modification  of  the  method  of  steepest  ascent. 

A  different  optimization  approach  based  on  the  concept  of  utility  or  desirability  was  followed 
by  Harrington  (1965)  and  Derringer  and  Suich  (1980).  In  this  approach,  each  response  function 
undergoes  a  certain  transformation  into  a  desirability  function,  <p,  such  that  0  s  ?  i  1.  The  choice 
of  the  transformation  depends  on  a  subjective  judgment  concerning  the  importance  (or  desirability) 
of  the  corresponding  response  values.  A  measure  of  the  overall  desirability  of  the  responses  is  ob¬ 
tained  by  combining  the  individual  desirability  functions  using  the  geometric  mean 
(<p,  x  ip3  x  •••  *  <t>,)'",  where  9,  is  the  desirability  function  for  the  rth  response.  The  simultaneous 
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optimization  of  the  r  response  functions  is  then  achieved  by  maximizing  the  overall  desirability 
function  over  the  experimental  region. 

More  recently,  Khuri  and  Conlon  (1981)  introduced  a  procedure  for  the  simultaneous  opti¬ 
mization  of  responses  that  are  represented  by  linear  models  of  the  form 

£(#)  *  .  f  ■  I,  2 . r , 

where  X  is  an  N  x  p  matrix  of  full-column  rank  and  &  is  a  vector  of  p  unknown  parameters 
(/'  *  1,  2,  ... ,  r).  The  rows  of  Y  ®  ...  :£,]  are  statistically  independent  with  each  having  a 

2ero  mean  and  a  variance-covariance  matrix  I.  A  distance  function  is  chosen  which  measures  the 
overall  closeness  of  the  response  functions  to  achieving  their  respective  optimal  values  at  the  same 
set  of  operating  conditions.  This  is  referred  to  as  an  'ideal  optimum*.  Optimum  operating  condi¬ 
tions  are  then  derived  by  maximizing  this  distance  function  over  the  experimental  region.  The 
distance  function  approach  permits  the  user  to  take  into  consideration  the  variance -covariance 
structure  for  the  r  responses. 

3.2  Sampling  Properties  of  Optima 

In  the  typical  analysis  of  a  second  degree  response  surface  one  is  interested  in 

i)  The  location  of  the  stationary  point 

ii)  The  response  at  the  stationary  point 

til)  The  characterization  of  the  stationary  point,  i.e.  as  a  point  of  maximum  or  minimum 
response  or  a  saddle  point 
iv)  The  contours  of  constant  response 

By  writing  the  true  response  equation  as 

ri  -  30  +  +  x!Bx  (3.4) 


where  s.  «  (jt„  Xj,  ... ,  x,)' ,  ft  -  (P„  p„ 
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The  response  at  the  stationary  point  is  given  by 

0.5) 

The  nature  of  the  stationary  point  can  be  ascertained  by  considering  the  matrix  B.  If  B  is  a  positive 
(negative)  definite  matrix  then  the  stationary  point  is  a  point  of  minimum  (maximum)  response. 
If  B  is  an  indefinite  matrix  the  stationary  point  is  a  saddle  point.  From  this  it  follows  that  the 
eigenvalues  of  B  are  important  to  the  proper  description  of  the  stationary  point. 

The  contours  of  constant  response  are  obtained  by  finding  all  values  of  ^  for  which  (3.4)  takes 
on  specified  values.  These  contours  are  useful  in  graphically  representing  a  response  surface  since 
they  permit  the  presentation  of  a  k  dimensional  surface  in  k  -  1  dimensions. 

Since  the  model  parameters  are  unknown,  they  must  be  estimated.  Point  estimates  of  the 
location  of  the  stationary  point,  the  response  there,  the  eigenvalues  of  B,  and  contours  of  constant 
response  can  be  obtained  by  replacing  the  model  parameters  by  their  respective  point  estimates 
obtained,  for  example,  by  the  methods  of  least  squares  or  maximum  likelihood.  Far  too  many  users 
of  RSM  allow  conclusions  to  be  drawn  concerning  the  nature  of  a  response  surface  and  the  location 
of  optimum  response  without  taking  into  account  the  distributional  properties  of  the  estimated  at¬ 
tributes  of  the  underlying  response  surface.  While  the  distribution  of  these  quantities  has  not  been 
considered  directly,  efforts  have  been  made  to  develop  interval  estimates  for  them. 
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the  location  of  the  stationary  point  can  be  shown  to  be 
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Box  and  Hunter  (1954)  consider  the  development  of  the  100(1  —  a)®/»  confidence  region  for 
the  location  of  the  stationary  point.  Letting  £  denote  the  coordinates  of  the  true  stationary  point 
and  assuming  the  adequacy  of  (3.4)  it  follows  that 

in. I-  a  Q  (3.6) 

d&  4 

If  5  represents  the  estimator  of  the  left-hand  side  of  (3.6)  formed  by  replacing  the  parameters  in  (3.5) 
with  normally  distributed,  unbiased  estimates,  it  follows  that 


where  the  elements  of  V  are  functions  of  the  elements  of  {XX)~ 1  and  X  is  associated  with  (3.4) 
when  expressed  in  general  linear  model  notation.  After  accounting  for  the  estimate  of  oJ  it  can  be 
shown  that 
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constitute  a  100(1  -  a)%  confidence  region  for  the  location  of  the  stationary  point.  Here  s1  is  the 
usual  error  mean  square.  It  should  be  noted  that  this  technique  is  applicable  to  any  surface  re¬ 
presented  by  a  linear  model  and  not  just  for  the  quadratic  model  considered. 

The  construction  of  a  confidence  interval  about  the  response  at  the  stationary  point  of  the 
true  surface  has  only  recently  been  a  subject  of  interest  in  the  statistical  literature.  Khuri  and 
Conlon  (1981)  give  an  expression  for  the  bounds  of  such  an  interval  conditional  on  the  estimated 
location  of  the  stationary  point.  For  the  general  linear  model 

l  *  i  ~  A'fQ.o2/) 


Spjotvoll  (1972)  and  Rao  (1973)  indicate  that  for  any  particular  continuous  function  gffi) 


Pr{min  $(fi)  *  *(&}  £  max  g(fi))  *  1  -  a 
gee  g  sc 
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where  C  -  {Q:  (fi  -  2)'*’*(2  -  Dips1  s  is  the  100(1  -  a)%  confidence  region  for 

0,  2  •  ( XX Y'X'z  is  the  least  squares  estimator  of  2.  *nd  J1  is  *■*  defined  previously.  Thus 
[mk  5(2),  m«  5(fi)]  is  a  confidence  interval  for  5(2)  with  a  confidence  coefficient  of  at  least 
(1  -  a).  Carter  et  al.  (1984)  make  use  of  this  result  to  obtain  a  conservative  confidence  interval 
on  the  response  at  the  stationary  point  in  RSM.  They  also  discuss  the  computational  aspects  of 
methodology. 

The  user  of  RSM  techniques  gains  considerable  insight  into  the  nature  of  the  stationary  point 
and  the  nature  of  the  response  system  from  estimates  of  the  eigenvalues  of  B.  As  a  result,  some 
indication  of  the  quality  of  estimation  of  the  eigenvalues  of  B  should  become  an  essential  part  of 
a  RSM  analysis.  Since  the  eigenvalues  of  the  B  matrix  can  be  expressed  as  particular  functions  of 
the  model  parameters,  the  approach  described  above  can  be  used  to  obtain  a  conservative 
100(1  -  a)%  confidence  interval  about  each  eigenvalue  of  B.  Carter  et  al.  (1986)  illustrate  this 
methodology  for  constructing  such  intervals  for  the  eigenvalues  of  B.  In  addition,  they  consider  the 
implications  of  having  such  an  interval  contain  zero  on  the  sampling  properties  of  the  estimated 
response  at  the  stationary  point.  They  conclude  in  this  case  that  a  change  in  the  strategy  of  analysis 
is  required  and  they  suggest  and  illustrate  one  such  change. 

It  is  sometimes  the  case  that  the  stationary  point  is  outside  the  region  of  interest  or  that  the 
stationary  point  is  unsatisfactory  due  to  associated  secondary  responses.  When  either  is  the  case, 
constraints  must  be  placed  on  the  optimization.  The  estimation  of  the  location  of  the  stationary 
point  subject  to  constraints  is  considered  elsewhere  in  this  review.  Stablcin  et  al.  (1983)  extend  the 
methodology  developed  by  Box  and  Hunter  (1954)  to  obtain  a  confidence  region  about  the  location 
of  the  stationary  point  in  the  presence  of  constraints.  In  their  ,  -  Carter  et  al.  (1986)  discuss  the 

construction  of  a  conservative  confidence  interval  about  the  response  at  a  constrained  stationary 
point.  It  is  shown  that  when  there  is  evidence  that  there  is  a  ridge  condition  as  opposed  to  a  true 
optimum  point,  i.e.,  the  confidence  interval  about  the  largest  or  smallest  eigenvalue  of  B  contains 
zero,  the  sampling  properties  of  the  response  at  an  unconstrained  optimum  may  be  undesirable. 
In  such  cases,  it  is  shown  that  use  of  a  region  constraint  as  in  a  ridge  analysis  will  provide  a  tighter 
confidence  interval  about  the  response  at  the  stationary  point. 


3.3  Use  of  Alternate  Models 


There  are  many  typical  RSM  situations  in  which  the  user  can  make  use  of  known  information 
concerning  the  distribution  of  the  observed  response.  Examples  of  such  response  variates  arc  sur¬ 
vival  time,  proportion  of  experimental  units  assigned  to  a  given  treatment  group  responding  in  a 
predetermined  way,  and  the  number  of  events  occurring  per  treatment  group,  i.e.  count  data.  To 
the  extent  that  a  continuous  symmetric  distribution  can  be  made  to  adequately  approximate  the 
distribution  of  such  data,  the  classical  application  of  RSM  has  provided  satisfactory  results.  Once 
a  proper  distributional  form  has  been  determined  the  usual  approach  is  to  express  the  mean  of  the 
distribution  as  a  function  of  the  experimental  conditions.  Most  often  a  second  order  polynomial 
is  used.  However,  in  their  review  article  Mead  and  Pike  (1975)  give  an  excellent  discussion  regarding 
the  choice  of  the  response  function  to  graduate  the  relationship  between  the  dependent  and  inde¬ 
pendent  variables.  They  provide  the  algebraic  form  of  a  number  of  different  response  functions  that 
have  been  used  and  discuss  the  usefulness  of  each.  The  choice  of  the  form  of  the  underlying  re¬ 
sponse  function  is  based  on  the  assumed  shape  of  the  underlying  relationship.  While  an  important 
topic,  it  is  limited  only  by  the  ability  of  the  investigator  to  formulate  the  mathematical  relationship 
with  the  desired  properties.  The  distributional  form  of  the  underlying  data  has  an  impact  on  the 
estimation  of  the  model  parameters  and  on  the  inferences  drawn  from  a  RSM  analysis.  Thus,  for 
the  remainder  of  this  section  some  developments  in  the  choice  of  a  distributional  form  will  be  re¬ 
viewed. 

Addelman  et  al.  (1966)  devised  an  experimental  scheme  for  estimating  the  optimal  combina¬ 
tion  of  two  drugs  in  preclinical  cancer  research.  The  procedure  is  such  that  at  most  three  exper¬ 
iments  must  be  performed.  The  first  experiment  uses  three  doses  from  within  the  effective  range 
of  each  of  the  individual  drugs  and  nine  combinations  arrived  at  by  considering  three  different  ratios 
between  the  levels  of  the  two  drugs.  The  data  from  this  experiment  are  analyzed  by  fitting  the 
quadratic  function 

y  ~  $ 0  +  +  ^2^2  7-  +  $22^2  + 

where 
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y  ■  survival  time  of  the  experimental  subject 

X,  ■  +  ^osc  ^cvcl  comP°un<*  *|  *  <  ■  1  •  2 

X,i  •  lower  extreme  of  effective  range  of  compound  1 

After  the  model  parameters  have  been  estimated  by  the  method  of  least  squares,  the  underlying 
response  surface  is  explored  and  the  optimal  combination  is  estimated.  From  this  a  second,  and 
sometimes  third,  experiment  is  performed  and  depending  upon  the  realization  of  given  conditions 
data  from  the  initial  experiments  are  included  in  the  final  analysis.  Unfortunately,  the  authors  il¬ 
lustrated  their  work  with  hypothetical  data.  Consequently,  it  was  not  possible  to  determine  the 
applicability  of  it. 

Due  to  the  nature  of  survival  data,  the  assumptions  required  for  the  appropriate  use  of  the 
method  of  least  squares  may  not  be  satisfied.  Cox  (1972)  proposed  the  proportional  hazards  model 
as  a  means  of  relating  survival  time  to  a  set  of  independent  variables.  Instead  of  working  directly 
with  survival  times  this  model  relates  the  hazard,  or  instantaneous  risk  of  failure,  k(r),  at  time  t  to 
the  concomitant  v  <  iables  z  as 

W)  -  Xo(/)exp  (*'B)  (3.7) 

where  Xq(()  is  the  hazard  function,  the  form  of  which  does  not  need  to  be  specified,  of  the  appro¬ 
priately  scaled  reference  group  and  2  is  the  vector  of  unknown  regression  coefficients  which  can  be 
estimated  by  the  method  of  maximum  partial  likelihood.  This  model  has  been  used  by  Carter  et 
al.  ( 1979)  in  preclinical  cancer  therapy  to  relate  the  risk  of  failure  to  dosage  levels  of  cytotoxic  agents 
given  in  combination.  Since  it  is  reasonable  to  assume  that  the  risk  of  failure  decreases  with  in¬ 
creasing  treatment  levels  until  a  minimum  is  reached  and  then  increases  due  to  toxicity  with  in¬ 
creasing  treatment  levels  beyond  the  minimum  point,  a  complete  second  order  polynomial  was  used 
to  approximate  in  equation  (3.7)  for  the  obvious  RSM  application.  An  exploration  of  the  fitted 
surface  followed. 

When  it  is  desired  to  relate  the  number  of  events  that  have  occurred  to  a  set  of  concomitant 
vanibles,  a  Foisson  regression  analysis  is  often  appropriate.  However,  in  many  experimental  situ¬ 
ations  counting  variables  exhibit  extra*  Poisson  variability  in  the  sense  that  there  is  more  variability 


in  the  data  set  than  can  be  explained  by  the  Poisson  distribution.  A  generalization  of  the  Poisson 
is  the  negative  binomial  distribution.  The  particular  form  of  the  negative  binomial  was  discussed 
by  Collings  and  Margolin  (1985).  Rac's  efficient  score  principle  (1973)  can  be  applied  to  test  the 
adequacy  of  the  Poisson  model.  Solana  et  al.  (in  press)  make  use  of  these  results  to  model  the  mean 
number  of  sister  chromatid  exchanges  induced  per  cell  as  a  polynomial  function  of  the  concen¬ 
trations  of  ethylnitrosourea  and  cis-diaminedichloroplatinum.  The  function  used  was  a  complete 
quadratic  with  three  additional  higher  order  interaction  terms.  Such  an  analysis  yielded  a  more 
complete  description  of  the  interaction  occurring  between  these  two  compounds  than  had  been 
obtained  previously. 

In  some  applications  binary  (e.g.  response/no  response)  responses  are  encountered.  In  such 
situations  it  is  often  desired  to  relate  the  proportion,  p,  of  responses  to  a  set  of  concomitant  vari¬ 
ables.  Assuming  a  constant  probability  of  response  within  treatment  groups,  this  has  beer  ac¬ 
complished  through  use  of  the  logit  transformation  of  the  binomial  parameter,  i.e., 

to  7~7  *  00  +  Mi  +  M2  +  -  +  0*** 

The  estimation  of  the  model  parameters  and  other  aspects  of  the  analysis  of  such  data  are  given 
by  Cox  (1977). 

Such  an  approach  has  been  used  to  relate  the  probability  of  a  favorable  treatment  outcome 
to  the  dosage  levels  of  two  anticancer  drugs  given  in  combination  (Carter  et  al.  1977).  The  use  of 
a  complete  second  order  model  was  justified.  The  fit  of  the  model  was  shown  to  be  adequate  and 
the  fitted  dose  response  surface  was  explored  and  the  estimated  optimal  treatment  was  obtained. 

It  sometimes  is  the  case  that  there  is  more  variability  in  the  data  than  can  be  explained  by  the 
binomial  distribution.  In  such  situations,  the  beta-binomial  distribution,  a  generalization  of  the 
binomial,  has  been  useful. 

The  form  of  the  beta  binomial  distribution  was  discussed  by  Segreti  and  Munson  (1981). 
As  in  the  use  of  the  negative  binomial  distribution  to  account  for  e;Ura-Poisson  variability  Rao's 
efficient  score  principle  (1973)  can  be  applied  to  test  the  adequacy  of  a  binomial  model.  Chinchilli 
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and  Chinnici  (1986)  have  used  this  approach  to  model  the  proportion  of  Drosophila  melanogasier 
eggs  that  hatch  as  a  function  of  exposure  to  caffeine  and  aflatoxin  Bl  given  in  combination.  They 
used  a  full  second  order  polynomial  to  approximate  the  dose-response  relationship  and  found 
extra-binomial  variability  in  their  data. 

IV.  SUBJECT  MATTER  APPLICATIONS 

In  the  two  previous  RSM  reviews,  much  attention  was  given  to  subject  matter  applications, 
v.  the  Hill  and  Hunter  review,  chemical  and  processing  applications  were  illustrated,  with  real  life 
examples  demonstrating  canonical  analysis  and  use  of  multiple  responses.  An  example  where  a 
canonical  analysis  provided  some  insight  ('elucidation')  into  the  mechanism  in  question  was  also 
given.  The  Mead- Pike  review  produced  examples  in  the  biological  area.  From  our  search,  it  is  clear 
that  standard  use  of  RSM  has  accelerated.  In  this  section  we  separate  applications  that  were  re¬ 
ported  in  subject  matter  journals  from  applications  that  are  finding  use  in  industry.  As  expected, 
many  imaginative  applications  were  found  in  the  chemical  and  processing  areas.  But  we  were 
somewhat  surprised  to  find  that  uses  have  spread  rather  dramatically  into  other  areas.  We  choose 
to  highlight  a  number  of  applications  here.  In  this  survey  of  the  subject  matter  literature,  several 
interesting  facts  become  apparent: 

i)  The  central  composite  is  used  more  than  any  other  family  of  RSM  designs. 

ii)  While  the  methods  of  exploring  the  response  surface  have  become  somewhat  more 
sophisticated,  users  nevertheless  are  not  close  to  making  maximal  use  of  the  enormous 
wealth  of  analytic  to  jIs  that  are  available.  This  is  due,  in  large  part,  to  the  *lag  time' 
that  naturally  exists  between  the  time  at  which  a  publication  appears  in  the  statistics 
literature  and  the  point  at  which  it  is  adopted  by  subject  matter  users.  However,  as 
in  the  case  of  many  areas  of  statistics,  our  communication  with  the  user  remains 
woefully  inadequate.  Perhaps  the  quickest  medium  of  communication  is  through  ef¬ 
fective  computer  software.  But  in  the  area  of  RSM.  there  has  been  little  planning  and 
leadership  in  this  area  of  development,  though  there  have  been  some  recent  rays  of 
hope. 
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iii)  Nearly  all  practical  RSM  problems  arc  truly  multiple  response  in  nature.  Many  users 
handle  multiple  responses  with  some  creativity  but  little  sophistication. 

In  what  follows,  examples  of  subject  matter  applications  are  given  in  areas  of  science  and  engi¬ 
neering,  biology  and  biomedicine,  and  applications  in  management  science,  operations  research,  and 
food  science. 

4. 1  Applications  in  Physical  and  Engineering  Sciences 

In  the  late  70's  scientists  dealing  in  the  important  area  of  gas  chromotography  began  using 
RSM  methods.  Scott  (1970),  Swingle  and  Rogers  (1972),  Kambara  et  al.  (1968),  and  Turina  et  al. 
(1974)  use  variations  of  RSM  to  determine  optimum  conditions  for  gas  chromotography.  Morgan 
and  Jacques  (1978)  study  relationships  between  flow  rate,  temperature,  stationary  phase  loading, 
and  the  responses,  level  of  separation  and  analysis  time.  They  indicate  how  helpful  an  under¬ 
standing  of  the  system  comes  from  elucidating  the  nature  of  the  response  surface.  Cross  sectional 
plots  in  this,  a  natural  three  variable  system,  were  used  to  estimate  the  conditions  where  desirable 
separation  and  analysis  time  were  achieved.  A  modified  simplex  algorithm  was  used  to  demonstrate 
the  minimization  of  analysis  time  at  fixed  threshold  values  of  peak  separation. 

Person  (1978)  used  RSM  notions  to  analyze  containment  leakage  from  a  sodium  Ore.  'Sen¬ 
sitivity  studies'  were  conducted  in  which  a  study  was  made  of  the  influence  of  sodium  mass,  con¬ 
tainment  volume,  radiation  heat  transfer  coefficient,  and  sodium  oxide  fraction  on  fraction  leakage. 
Prato  and  Moms  (1984)  studied  the  effect  of  detergent  concentration,  agitation  time,  and  washing 
time  on  the  amount  of  soil  removed  from  fabric.  In  order  to  underscore  the  natural  application 
of  RSM  and  second  order  composite  designs  in  human  factors  engineering  type  experiments,  Simon 
(1970)  illustrates  with  a  study  involving  three  display  panel-type  variables  and  their  effect  on  target 
recognition  for  airline  pilots. 

Oiivi  (1980)  discusses  the  need  for  use  of  RSM  in  exploring  and  identifying  certain  features 
of  systems.  He  uses  as  a  'typical'  application  a  study  involving  five  design  variables  with  an 
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orthogonal  central  composite  design.  The  application  involves  factors  that  influence  ballooning 
time,  an  important  variable  in  nuclear  safety. 

Bodden  and  Edwards  (1982)  use  a  Box-Bchnken  design  in  an  RSM  study  in  the  investigation 
of  the  mechanism  involved  in  the  assay  of  cretinine.  Claycomb  and  Sullivan  (1976)  use  a  three 
factor  central  composite  design  and  a  ridge  analysis  of  the  data  to  illustrate  the  methodology  for 
selecting  a  cutting  tool  for  maximization  of  profit. 

Contour  plots  of  constant  response  without  an  analytic  method  for  finding  optimum  condi¬ 
tions  are  often  the  source  of  conclusions  drawn  by  RSM  analysts.  Such  was  the  case  in  the  work 
by  Bretaudiere  et  al.  (1980)  in  a  study  of  the  Analine  Aminotransferase  Kinetic  Assay.  Fujiwara 
et  al.  (1975)  study  the  atomization  mechanism  in  air-acetylene  flames  using  response  contours. 
Burtis  et  al.  (1981)  used  data  from  a  central  composite  design  to  produce  useful  contours  of  con¬ 
stant  response  and  a  canonical  analysis  that  allowed  for  optimization  of  a  kinetic  method  that  ap¬ 
plies  to  the  enzymatic  measurement  of  ethanol.  Response  contours  are  also  the  subject  of 
information  obtained  from  a  response  surface  analysis  in  which  interest  centers  on  formulation  and 
optimization  of  vitreous  bonds  in  a  paper  by  Chait  and  Fotlz  (1981).  The  field  of  nuclear  engi¬ 
neering  has  not  been  without  studies  involving  RSM,  particularly  where  simulation  is  involved. 
Heller  et  al.  (1977)  use  data  from  a  central  composite  design  to  study  the  mechanism  involved  in 
thermal-hydraulic  margin  analysis.  Vaurio  and  Mueller  (1978)  use  response  surfaces  generated  from 
reactor  simulations  to  obtain  probability  distributions  of  selected  consequences  of  a  liquid-metal 
fast  breeder  reactor  core  accident. 

4.2  Applications  in  Food  Science 

The  food  industry  has  been  a  prime  user  of  RSM  since  the  early  1970's.  This  is  reflected  in 
study  of  content  of  journal  articles  and  industrial  usage  as  well.  We  cannot  begin  to  cite  ail  of  the 
work  in  this  area.  Johnson  and  Zabik  (1981)  use  a  mixture  design  to  build  a  response  surface  to 
study  interactions  among  proteins  in  angel  food  cake.  Lah  et  al.  (1980)  use  a  27'J  factorial  with 
an  impressive  application  of  steepest  ascent  to  optiro.,.  whipping  properties  of  an  ultraflltered  soy 
product.  R.G.  Henika  began  using  RSM  at  Foremost  Foods  in  the  early  1970's.  His  leadership 
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resulted  in  usage  by  others  even  though  many  of  the  applications  in  this  field  involve  rather  pre¬ 
carious  or  tricky  panel  data  from  sensory  evaluations.  The  papers  by  Henika  and  Palmer  (1976) 
and  Henika  (1982)  are  certainly  important  contributions.  In  the  latter,  particular  attention  is  paid 
to  applications  with  sensory  data. 

Second  order  models  were  used  to  find  the  nature  of  the  system  relating  gelling  properties  of 
a  meat  /oa/ analog  to  pH,  salt  content,  cooking  time,  soil  isolate,  and  egg  white  teve?.  While  formal 
optimization  procedures  were  not  used,  the  response  surface  allowed  Jao  et  al.  (1982)  to  better 
understand  the  mechanism.  Desirable  gelling  properties  can  be  obtained  with  reduced  egg  white 
usage  with  proper  manipulation  of  salt  content  and  cooking  time. 

The  extrusion  of  starch  has  important  applications  in  food  and  industrial  applications.  El* 
Dash  et  al.  (1983)  used  response  surface  methods  to  find  ideal  operating  conditions  in  this  area  of 
technology.  They  discovered  that,  depending  on  the  application,  one  can  control  levels  of  extrusion 
temperature,  screw  speed,  and  moisture  content  to  produce  the  proper  level  of  gelatination. 

Daley  et  al.  (1978)  used  RSM  to  determine  the  optimal  combination  of  soy,  sodium 
tripolyphosphate,  and  water  to  produce  the  most  acceptable  sausage-type  product  from  minced 
mullet.  A  central  composite  design  was  used  to  study  the  effects  of  calcium  chloride  and  cysteine 
on  heat  induced  whey  protein  gelatin  by  Schmidt  et  al.  (1979).  After  the  parameters  of  a  second 
oxdeT  polynomial  were  estimated,  plots  of  the  fitted  surface  were  generated.  Sefa-Dedeh  and  Stanley 
(1979)  made  use  of  a  central  composite  design  to  fit  the  percent  of  nitrogen  extracted  from  a  type 
of  meal  as  a  complete  quadratic  in  time  of  extraction,  temperature  of  extraction  and  meal  to  solvent 
ratio.  The  fitted  polynomial  was  explored  to  estimate  the  levels  of  the  independent  variables  asso¬ 
ciated  with  optimal  extraction  of  protein  nitrogen.  Mixture  response  surface  methodology  was  used 
by  Huor  et  al.  (1980)  in  conjunction  with  laboratory  sensory  evduation  and  small  scale  consumer 
tests  to  optimize  the  proportions  of  watermelon,  pineapple,  and  orange  juice  in  a  fruit  punch.  Min 
and  Thomas  (1980)  utilized  RSM  analytic  techniques  to  determine  the  relationship  between  ingre¬ 
dients  and  physical  characteristics  of  dairy  whipped  topping  and  optimize  the  ingredient  concen¬ 
tration  by  the  simultaneous  analysis  of  fat,  com  syrup  solids,  and  stabilizer  to  produce  a  dairy 
hipped  topping  which  could  be  frozen  for  storage  and  distribution  without  adverse  affects.  Lee 
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and  Hoscncy  (1982)  used  response  surface  methodology  to  optimize  the  formulation  of  single-stage 
cake  mixes  for  White  layer  cakes.  RSM  techniques  were  used  by  Tseo  et  al.  (1983)  to  determine 
the  optimum  combination  of  levels  of  washing  temperature,  washing  ratio  of  water  volume  to 
sample  weight  and  washing  time  on  the  quality  of  minced  mullet  flesh.  Tong  et  al.  (1984)  used  a 
central  composite  design  to  study  the  effects  of  safflower  oil  concentration,  emulsifier  concentration 
and  freezing  temperature  on  maximum  overrun  and  fat  destabilization  of  ice  cream.  McLellan  et 
al.  (1984)  used  a  rotatable  experimental  design  to  generate  data  for  a  sensory  analysis  of  carbonated 
apple  juice  as  a  function  of  levels  of  carbonation  and  soluble  solids.  The  fitted  surface  for  accepta¬ 
bility  was  described  as  having  a  ridge  of  high  acceptability  running  through  the  range  of  typical 
soluble  solids  and  carbonation  levels. 

4.3  Applications  in  Social  Sciences 

Much  use  has  been  made  of  RSM  in  certain  areas  in  social  sciences.  Economics,  operations 
research,  and  system  simulation  are  but  a  few  fields  that  have  benefitted.  Shechtcr  and  Heady 
(1970)  use  response  surface  techniques  to  design  and  analyze  experiments  from  a  simulation  model 
dealing  with  the  feed  grain  program.  Four  responses,  net  farm  revenue,  net  farm  revenue  partic¬ 
ipants,  stock  accumulation,  and  government  costs  were  considered.  Multiple  response  analysis  is 
discussed.  In  this  case,  the  analysis  allows  for  an  optimum  decision  rule  regarding  trade-offs  be¬ 
tween  increasing  farm  income  and  reducing  government  costs.  Montgomery  and  Bettencourt 
(1977)  produced  an  excellent  account  in  which  the  advantages  of  RSM  methods  in  simulation 
studies  were  reviewed.  Heavy  emphasis  is  put  on  use  in  multiple  response  studies.  An  interesting 
illustration  was  given  in  which  a  simulation  of  a  military  tank  duel  is  analyzed  to  ascertain  the 
values  of  two  design  variables  (mean  time  to  fire  first  round  and  mean  time  between  rounds)  that 
give  desirable  values  of  four  responses.  They  use  a  nonlinear  programming  technique  to  analyze 
data  taken  from  a  rotatable  central  composite  design. 

Montgomery  and  Evans  (1975)  discuss  the  use  of  various  classes  of  second  order  response 
surface  designs  in  simulation  work.  Two  variable  equiradial  designs,  uniform  precision  rotatable 
ccd  s  and  other  designs  are  discussed.  An  illustration  is  given  in  which  a  simulation  of  a  six  inter- 


section  street  network  system  was  used.  Data  was  taken  using  a  rotatable  c.c.d.  with  two  design 
variables.  Canonical  analysis  revealed  a  maximum  at  the  stationary  point.  The  results  produced 
optimum  conditions  on  the  selection  of  traffic  signal  settings.  Hunter  and  Naylor  (1970)  point  out 
that  a  simulation  system  is,  indeed,  the  source  of  an  experiment  similar  to  a  biological  or  chemical 
system  and  thus  experimental  design  and  RSM  methods  do  apply.  Factorial  and  fractional  factorial 
arrays  are  discussed  in  the  context  of  an  inventory  problem. 

Smith  (1975)  made  an  empirical  study  of  various  optimum  seeking  procedures  that  are  fine* 
quently  used  in  computer  simulation  situations.  The  methods  studied  were  random  search,  one 
factor  at  a  time  movement,  RSM  Version  /  and  RSM  Version  II.  Version  1  involves  use  of  two- 
level  first  order  designs  and  the  method  of  steepest  ascent,  while  Version  II  involved  steepest  ascent 
with  the  use  of  the  simplex  design.  In  addition,  acceleration  versions  of  steepest  ascent  were 
studied.  RSM  Version  I  appeared  to  give  best  results.  Any  attempt  to  accelerate  the  steepest  ascent 
procedure  was  fruitless.  Biles  (1981)  describes  a  two  phase  procedure  in  which  a  complex  search 
procedure  is  employed  followed  by  a  second  order  analysis  of  the  resulting  data.  An  illustration 
of  the  method  is  given  from  an  inventory  system. 

4.4  Applications  in  the  Biological  Sciences 

The  techniques  that  comprise  a  response  surface  analysis  are  being  used  more  and  more  often 
in  the  biological  sciences.  This  is  evidenced  by  the  large  number  of  citations  in  the  bibliography 
of  this  report.  In  this  section  some  of  the  more  novel  applications  will  be  highlighted. 

RSM  techniques  have  been  found  useful  in  the  study  of  the  relationship  between  the  chemical 
structure  of  a  compound  and  its  biological  activity.  Mager  (1982a,  1983)  studied  the  structure- 
neurotoxicity  relationship  of  organophosphorous  pesticides  in  this  manner  and  used  a  canonical 
analysis  of  the  fitted  equation  to  elucidate  properties  of  the  response  surface.  Mager  (1982b)  also 
studied  the  activity  of  rifamycin  derivatives  against  a  rifaraycin-MS-resistant  mutant  of 
Staphylococcus  aureus. 

Several  references  in  the  microbiologic  literature,  Maddox  and  Richert  (1977),  Cinto  et  al. 
(1977)  and  Farrand  et  al.  (1983),  indicate  the  useful  application  of  RSM  techniques,  including  the 
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use  of  rotatable  central  composite  designs,  to  optimize  microbiological  media  for  the  growth  of 
various  microorganisms.  Cheynier  et  aj.  (1983)  were  able  to  demonstrate  that  a  yeast  isolated  from 
the  digestive  tube  of  the  larva  of  a  parasite  of  eucalyptus  trees  was  capable  of  bioconverting 
citroncllal  to  citronellol.  RSM  procedures  were  used  to  achieve  the  optimization  of  the  exper¬ 
imental  conditions  for  that  byconversion  process. 

Various  industrial  pollution  studies  have  employed  response  surface  methodology.  Huck  et 
al.  (1977)  determined  the  polymer  properties  and  mixing  conditions  required  to  produce  optimal 
flocculation  for  mine  waters  of  specified  strengths  containing  iron,  zinc,  and  copper  either  singly 
or  in  combination.  Wallis  (1978)  reports  on  the  use  of  RSM  in  studies  related  directly  to  power 
station  cooling  systems. 

Dincer  and  Ozdurmus  (1977)  used  the  method  of  steepest  descent  to  determine  the  most 
suitable  combination  of  four  independent  formulation  and  process  variables  for  the  disintegration 
time  of  coated  tablets  in  simulated  intestinal  fluid.  Shek  et  al.  (1980)  evaluated  the  potential  of  the 
Nelder-Mead  simplex  search  procedure  for  optimizing  a  capsule  formulation.  At  the  completion 
of  the  search,  these  investigators  fitted  a  polynomial  model  to  the  data  and  plotted  the  estimated 
response  surface.  In  a  similar  study  Chowhan  et  al.  (1982)  studied  the  effects  of  moisture  and 
crushing  strength  on  tablet  friability  and  in  vitro  dissolution.  They  used  a  complete  second  order 
model  in  the  independent  variables  to  graduate  the  response  surface.  The  fitted  model  was  explored 
by  plotting  contours  of  constant  response  and  estimating  the  location  of  the  stationary  point.  Fast 
et  al.  (1983)  used  a  central  composite  design  with  five  independent  variables  in  an  assay  for  creatine 
kinase.  A  simplex  maximization  algorithm  was  used  to  determine  areas  of  maximum  sensitivity. 
Bclloto  et  al.  (1985)  used  RSM  to  study  the  solubility  of  pharmaceutical  formulations.  They  also 
discussed  the  optimization  of  one  response  in  the  presence  of  constraints  due  to  the  consideration 
of  other  response  variables. 

There  are  numerous  reports  in  the  clinical  chemistry  literature  on  the  use  of  response  surface 
methods,  Rautela  et  al.  (1979)  describe  the  theory  and  application  of  the  response  surface  approach 
to  simultaneous  optimization  of  multiple  interdependent  variables.  .They  point  out  that  such  an 
approach  permits  the  determination  of  accurate  optima  which  is  required  for  the  formulation  of 


58 


analytically  reliable  clinical  methods.  Thompson  et  al.  (1981)  report  on  the  use  of  simultaneous 
RSM  optimization  techniques  to  optimize  an  assay  for  alkaline  phosphatase.  London  et  al.  (1982)  - 
use  RSM  to  optimize  the  assay  of  gamma-glutamyltransferase.  The  authors  note  that  while  such 
techniques  do  not  yield  a  mechanistic  understanding  of  an  enzyme  assay  they  do  produce  an  oper- 
ational  understanding  of  how  an  assay  functions.  In  a  similar  paper,  Coleman  et  al.  (1983)  use 
RSM  to  optimize  assays  for  antithrombin  III  and  plasminogen. 

Roush  et  al.  (1979)  discuss  the  usefulness  of  a  RSM  approach  in  poultry  nutrition  research. 
Body  weight  and  feed  conversion  response  for  Japanese  quail  were  optimized  as  functions  of  protein 
and  energy  levels.  Heady  et  al.  (1980)  express  egg  production  as  a  quadratic  function  of  the 
amounts  of  com  and  soybean  meal  with  time  as  an  independent  variable.  Contours  of  constant 
response  were  generated  to  help  determine  the  least  cost  feed  mix.  Roush  (1983)  utilizes  a  central 
composite  experimental  design  with  male  broilers  to  examine  quantitatively  the  protein  levels  in 
starter  and  finisher  rations  and  the  time  of  ration  change  to  optimize  body  weight,  carcass  weight, 
feed  conversion,  and  net  profit.  The  optimization  was  done  for  each  response  variate  independently 
of  the  other  responses. 

Response  surface  methods  have  been  used  to  elucidate  the  actions  and  interactions  of 
cytotoxic  drugs  in  combination  and  to  estimate  the  optimal  levels  of  each  drug  for  the  treatment 
of  cancer  with  and  without  side  effect  constraints.  Carter  et  al.  (1983)  discuss  the  use  of  logistic  and 
proportional  hazards  regression  models  to  this  end.  These  authors  give  software  for  the  estimation 
of  the  stationary  point  and  confidence  regions  about  its  location  obtained  from  the  use  of  these 
models.  The  history  of  the  use  of  RSM  in  cancer  research  has  been  reviewed  by  Carter  and 
Wampler  (1986).  While  the  examples  used  involve  data  from  animal  studies  there  is  interest  in 
applying  these  techniques  to  clinical  studies.  The  biggest  change  required  to  accomplish  this  will 
be  in  the  design  of  the  studies.  The  classical  clinical  trial  in  cancer  research  results  in  a  number  of 
patients  being  randomized  to  two  or  more  treatment  groups  and  followed  until  an  event  occurs. 
The  number  of  patients  assigned  is  determined  so  that  there  is  a  reasonable  chance  of  detecting  a 
difference  of  a  predetermined  size  among  the  treatment  groups.  While  such  an  approach  provides 
considerable  information  about  a  few  treatment  groups,  it  does  not  offer  much  about  the  underlying 


dose  response  relationship.  Box  (1958)  in  a  discussion  of  a  paper  by  Mantel  (1958)  was  early  to 
suggest  the  usefulness  of  RSM  in  this  area.  Currently  the  Mid  Atlantic  Oncology  Program  is  run¬ 
ning  a  clinical  trial  in  the  study  of  breast  cancer  that  is  designed  as  a  factorial  experiment  and  will 
be  analyzed  using  response  surface  techniques. 

4.5  Industrial  Use  of  RSM 

The  1980's  brought  a  new  urgency  in  industry  in  the  U.S.,  an  urgency  that  underscores  the 
need  for  quality.  There  is  a  new  motivation  in  American  industry  brought  about  by  the  success 
enjoyed  in  Japan  in  the  use  of  Quality  Control  and  efficient  statistical  methods.  The  areas  of  ex¬ 
perimental  design  and,  specifically,  response  surface  methods  are  receiving  considerably  more  at¬ 
tention  than  ever  before.  This  new  and  renewed  interest  are  'fallout'  from  this  search  for  quality 
and  precision. 

In  our  communication  with  industrial  statisticians  and  research  workers,  we  found  use  of 
RSM  among  a  wide  variety  of  types  of  industrial  companies.  In  most  cases  they  sought  and,  of 
course,  were  promised  anonymity.  As  expected,  we  found  response  surface  analysis  and  design  in 
the  chemical  industry,  and  we  were  able  to  uncover  a  growing  number  of  applications  in  the  areas 
of  foods,  tobaccos,  military  research,  pharmaceuticals,  petroleum,  electronics,  and  many  other 
fields.  Our  search  did  reveal  some  interesting  but  not  unexpected  information: 

i)  Aside  from  the  central  composite  and  uox-Behnken  designs,  fractions  of  2-level  fac¬ 
torials  and  Plackett-Burman  designs  are  used  to  a  large  extent,  with  the  simplex  design 
used  in  conjunction  with  the  Nelder-Mead  Simplex  Search  Algorithm.  There  is  also 
a  growing  usage  of  computer  generated  RSM  designs. 

ii)  In  regard  to  RSM  analysis,  some  analysts  cany  the  computations  through  to  an  opti¬ 
mization  phase,  while  others  do  not.  Ridge  analysis  is  used  for  finding  candidates  for 
optimum  operating  conditions  though  usage  is  not  as  extensive  as  one  might  expect. 
Some  users  merely  use  graphical  methods  with  overlaid  contour  plots. 


There  is  no  universal  software  tool  for  doing  RSM  analysis.  In  particular,  some  users  use 
1MSL  (International  Mathematical  and  Statistical  Libraries)  subroutines,  while  some  resort  to  in* 
house  software  versions  for  finding  optimum  conditions.  Some  users  in  industry  appear  frustrated 
by  the  lack  of  software  for  finding  optimum  conditions  and  resort  to  two  dimensional  contours  for 
analysis. 

In  what  follows,  a  sample  of  RSM  scenarios  is  given  that  apply  to  actual  RSM  users  in  in* 
dustry.  We  have  attempted  here  to  give  samples  across  many  fields  of  application. 

A  tobacco  company  uses  RSM  as  a  primary  tool  for  elucidation  of  fundamental  relationships 
and  searching  for  optimum  conditions  in  the  case  of  nearly  every  new  product  considered  in  a  re¬ 
search  effort.  They  use  no  software  for  design  construction  but  use  in-house  software  developed 
from  numerical  analysis  algorithms  for  analysis.  The  latter  allows  for  optimization  under  various 
types  of  constraints  on  the  design  variables.  The  designs  used  vary  considerably.  For  problems 
where  cost  is  important,  Plackett-Bunnan  screening  designs  are  used  and  hybrid  and  small  com¬ 
posite  designs  are  used  in  the  second  order  case.  Central  composite  and  Box-Behnken  designs  are 
used  extensively.  Problems  are  almost  always  multiple  response  in  nature,  though  the  responses 
are  too  many  in  number  to  often  do  formal  multiple  response  optimization.  Some  of  the  responses 
are  discrete  'scores'  from  sensory  analysis  but  many  responses  are  naturally  continuous  in  nature. 

A  large  oil  company  uses  RSM  to  a  moderate  extent.  Applications  vary  from  polymer  opti¬ 
mization  to  the  exploration  of  a  detergent  system.  In  some  applications  the  design  variables  are 
of  the  mixture  type.  Mixture  variables  and  process  variables  often  occur  in  the  same  setting.  The 
computer  program,  ACED  (see  Welch  1985),  is  often  used  to  generate  the  design  in  the  mixture 
situation.  Basic  2-level  factorials  and  central  composite  designs  are  also  used.  In  some  cases,  de¬ 
termination  of  optimum  conditions  is  important.  The  algorithms  in  GRG2  (see  Cornell  1981)  and 
XSTAT  (1984)  are  generally  used  to  find  optimum  conditions.  This  company  also  envisions  con¬ 
sidering  RSM  in  the  future  for  plant  scale  problems  dealing  in  perhaps  more  than  1000  variables, 
i.e.,  allow  the  technique  to  play  a  role  in  control  theory. 

A  chemical  company  finds  that  RSM  design  and  analysis  can  be  used  to  solve  two  types  of 
problems.  The  first  is  a  region  seeking  procedure  in  which  one  is  seeking  the  general  vicinity  of 


•best'  operating  conditions.  The  method  of  EVOP  described  by  Box  (1957)  is  used  with  slight 
modification.  Steepest  ascent  has  been  used  at  times  with  little  success.  Some  success  has  been 
experienced  with  the  Spendley  et  al.  (1962)  simplex  algorithm.  A  computer  program  has  been 
developed  which  makes  use  of  various  types  of  designs,  allows  input  of  experimental  results,  and 
calculation  of  the  next  experimental  run.  It  has  been  found  that  the  most  difficult  aspect  of  this 
procedure  is  obtaining  from  the  experimenter  a  commitment  on  combating  response  values  from 
the  multiple  dependent  variables  so  that  the  simplex  will  move.  In  the  second  problem  encountered 
by  the  chemical  company,  they  are  near  the  optimum  and  wish  to  locate  or  confirm  its  exact  lo¬ 
cation.  Central  composite  designs  are  used  almost  exclusively.  Small  composite  designs  are  often 
used  with  no  concern  for  cither  orthogonality  or  rotatability.  An  in-house  program  is  used  to 
construct  central  composite  designs  for  researchers.  The  model  fitting  aspect  of  the  analysis  is  done 
by  the  all  possible  subset  procedure  in  BMDP.  Optimization  is  'accomplished'  with  two  dimen¬ 
sional  contour  plots. 

A  large  food  company  uses  RSM  in  several  cases.  It  is  used  in  product  development  in  de¬ 
veloping  new  products  and  in  cost  reducing  old  products.  The  responses  are  either  continuous 
variables  or  consumer  acceptability  ratings.  RSM  is  also  used  in  process  control  and  in  the  evalu¬ 
ation  of  new  equipment.  Generally,  for  first  order  models,  2-level  factorial  or  Plackett-Bunnan 
designs  are  used.  In  the  case  of  second  order  models,  Box-Behnken  designs  and  central  composite 
designs  are  used.  At  times,  the  c.c.d.  is  resisted  by  the  experimenter  because  five  levels  become  in¬ 
convenient  because  it  is  difficult  to  adjust  equipment  accordingly.  Errors  in  controlling  factor  levels 
occur  quite  often. 

Another  large  and  well  diversified  chemical  company  makes  use  of  RSM  on  a  research  level. 
This  company  has  a  large  group  of  professional  statisticians  and  an  in-depth  program  to  teach  basic 
experimental  design  to  engineers  and  scientists.  They  make  use  of  an  in-house  computer  program 
that  assists  non-statisticians  in  designing  and  analyzing  experiments.  The  program  involves  contour 
plots  and  optimization  routines.  There  is  great  potential  for  use  of  RSM  in  this  company,  and 
most  problems  are  multiresponsc  in  nature.  Mixture  problems  are  the  rule,  and  an  algorithm  has 
been  programmed  in-house  to  generate  mixture  designs  described  in  Thompson  and  Myers  (1968). 
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Extreme  vertices  designs  are  used  as  well.  In  addition,  considerable  use  is  made  of  the  ACED 
(Welch,  1985)  algorithm  for  design  construction. 

A  large  research  and  development  organization  makes  use  of  RSM  in  several  topical  areas, 
including  chemical  and  fermentation  optimization  and  development  of  consumer  food  products. 
In  addition,  it  is  also  used  for  optimization  of  settings  associated  with  processes  in  the  refining  and 
hydrogenation  of  vegetable  oils.  In  almost  all  cases,  'optimization*  requires  amving  at  the  most 
agreeable  compromise  among  various  conflicting  responses.  In  the  case  of  two  design  variables,  the 
choice  of  best  conditions  is  found  by  use  of  graphical  overlay  contour  plots.  For  more  than  two 
variables,  a  utility  function  is  constructed  and  a  Nclder-Mead  simplex  optimization  routine  is  used 
to  identify  an  acceptable  region  in  the  design  space,  in  some  instances  ridge  analysis  is  used  to  ei* 
ther  find  optimum  conditions  in  a  single  response  situation  or  to  find  regions  for  subsequent  ex¬ 
perimentation.  The  class  of  central  composite  designs  is  relied  upon  heavily.  Blocking  is  often 
required  and  c.c.d.'s  are  used  which  block  orthogonally.  At  times  rather  involved  constraints  are 
induced  on  the  'optimum'  conditions  la  this  case  they  use  'brute  force'  grid  search  routines  that 
are  supplemented  with  confirmatory  inais. 

A  consulting  firm  used  the  results  of  an  RSM  study  to  develop  a  generic  method  for  evalu¬ 
ating  submarine  tracking  algorithms.  Simulation  experiments  were  conducted  under  common  sets 
of  scenarios  and  then  algorithm  performance  was  evaluated  according  to  the  behavior  of  the  re¬ 
sponse  surfaces  generated.  RSM  is  also  under  consideration  for  several  sonar  operator  studies  to 
determine  the  effects  on  performance  of  various  processing  schemes.  Their  analyses  are  performed 
using  IMSL  subroutines. 

V.  FUTURE  DIRECTIONS  OF  RSM 

To  comment  about  future  directions  of  RSM,  one  must  ponder  the  present  status  with  regard 
to  utilization  and  importance.  It  is  rather  unfortunate  that  the  expanse  of  knowledge  gained  over 
the  past  25  years  in  RSM  has  not  blossomed  into  a  full-scale  deployment  of  the  techniques  devel¬ 
oped  in  this  area.  As  we  indicated  earlier,  RSM  is  being  used  in  industry  more  and  more,  but  the 
use  does  not  reflect  recent  advances  as  much  as  one  might  wish.  This  is  true  more  in  analysis  than 
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design,  The  'fallout'  of  the  recent,  almost  passionate,  'push'  for  quality  and  experimental  design 
in  industry  is  bringing  RSM  to  the  attention  of  many  potential  users.  However,  usage  lags  far  be* 
hind  its  potential.  Several  factors  may  have  contributed  to  this,  the  most  prominent  of  which  we 
believe  are  explained  in  the  following  paragraphs. 

There  is  a  severe  lack  of  communication  between  academic  statisticians  in  this  country  and 
those  in  industry.  There  are  some  attempts  at  communication  and  few  successes.  This  is  true  in 
many  areas,  not  merely  RSM.  One  testament  to  this  is  the  fact  that  Japanese  statisticians  give 
lectures  and  short  courses  to  industrial  companies  in  this  country  on  effective  use  of  experimental 
designs  that  are  based,  in  large  part,  on  basic  concepts  developed  in  the  U.S.  yean  ago.  There  is 
no  question  that  the  historical  underpinnings  of  RSM  were  produced  with  particular  applications 
in  mind  When  we  move  away  from  reseat  ch  that  does  not  have  particular  problem  solving  as 
motivation,  we  lose  our  line  of  communication.  Statisticians  in  academia  bear  much  but  not  the 
entire  burden.  They  have  little  access  to  what  is  actually  going  on  in  industry  as  far  as  the  de¬ 
ployment  (or  lack  of  it)  of  RSM  in  concerned.  This  can  be  partly  attributed  to  the  fact  that 
'sometimes,  those  in  industry  are  dissuaded  from  publishing  for  fear  of  revealing  proprietary  infor¬ 
mation',  as  Gerald  Hahn  (1984)  explains.  As  a  result,  some  of  the  research  that  has  been  under¬ 
taken  in  RSM  was  not  geared  to  solving  a  particular  practical  problem.  Even  a  search  of  the 
literature,  statistical  or  otherwise,  in  an  effort  to  find  a  real  data  set  to  use  in  an  example,  can 
sometimes  be  unwieldy. 

There  is  a  severe  software  problem.  It  is  unfortunate  but  true  that  data  analysts  are  attracted 
to  statistical  techni  ques  for  which  software  has  been  written.  There  are  many  historical  illustrations 
of  this  in  our  field.  Good  efficient,  well  documented  software  is  a  swift  and  sure  form  of  commu¬ 
nication.  As  we  indicated  in  a  previous  section,  there  is  no  commercial  up-to-date  RSM  package 
for  doing  analysis.  Of  course,  there  are  many  new  avafiabie  packages  for  constructing  experimental 
designs,  and  RSM  designs  are  the  subjects  of  some  of  these  packages.  But  as  it  was  expressed  in 
the  Steinberg  and  Hunter  experimental  design  review,  the  'expert  system'  for  design  generation  is 
not  near,  particularly  in  RSM  situations. 
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There  is  apathy  on  the  part  of  some  practitioners  who  either  believe  that  RSM  is  of  little  use 
to  them,  or  simply  ignore  new  developments  in  favor  of  antiquated,  and  often  ad  hoc,  methods. 
There  is  also  insufficient  training  in  RSM.  In  academia,  very  few  statistics  departments  offer  a 
course  in  RSM.  In  industry,  little  exposure,  if  any,  to  RSM  is  provided  in  the  form  of  short 
courses. 

The  review  paper  by  Steinberg  and  Hunter  (1984)  along  with  the  several  discussion  papers 
that  followed  it  have  expressed  concern  about  the  above-mentioned  communication  problems  and 
offered  several  remedies.  We  believe  that  greater  emphasis  should  be  placed  on  the  teaching  of 
RSM.  A  belter  cooperation  and  more  formal  lines  of  communication  between  academia  and  in¬ 
dustry  is  needed.  Statisticians,  scientists,  and  engineers  can  learn  and  benefit  from  one  another. 
This  fact  is  certainly  not  specific  to  RSM  but  applies  more  broadly  in  our  profession.  Far  too  many 
academicians  conduct  and  are  rewarded  (i.e.,  tenure  ard  promotion)  for  research  in  which  math- 
errat'  *legance  is  far  more  prominent  than  application  to  real  problems.  More  academic  re¬ 
st  need  to  let  their  work  be  motivated  by  real  problems,  much  like  the  case  of  Box  in  the 

early  50's  and  Scheffe  with  the  development  of  mixture  designs  in  the  mid  50's.  New  developments 
in  RSM  should  be  made  accessible  to  practitioners  and  research  workers  with  the  proper  software. 
In  return,  engineers  and  applied  statisticians  working  in  industry  should  make  a  genuine  clTort  to 
share  their  experience  with  those  in  academia. 

There  is  no  doubt  that  many  scientists  and  engineers  can  benefit  from  tools  of  RSM  but  do 
not  know  that  they  exist^!  We  do  not  presume  that  this  is  a  problem  solely  with  RSM.  Our  earlier 
comments  may  very  well  be  relevant  in  regard  to  ail  areas  of  applied  statistics.  However,  the 
enormity  of  the  problem  becomes  evident  when  one  'reaches  out'  and  essentially  contacts  and 
interviews  many  scientists  in  industry  as  we  have  done.  In  what  follows,  we  indicate  what,  in  our 
opinion,  are  proper  directions  for  future  research.  But  let  us  emphasize  that  users  are  sufficiently 
"far  behind'  in  use  of  RSM  tools  that  the  need  for  proper  communication  is  far  more  pressing. 


Future  Research  Directions 

Multiple  responses,  sequential  designs  and  analyses,  and  designs  for  the  less  orthodox  RSM 
applications  are  many  areas  where  one  may  benefit  from  immediate  research.  Both  design  and 
analysis  involving  multiple  responses  have  received  precious  little  attention  in  spite  of  the  acute  need 
for  it.  Biologists  as  well  as  engineers  are  finding  applications  of  RSM  in  which  the  model  is  non* 
linear,  e.g.,  logistic  and  other  standard  growth  models.  A  fruitful  area  of  research  deals  in  the  study 
of  designs  or  certain  specific  nonlinear  functional  forms. 

Much  has  been  written,  of  course,  about  computer  generated  design.  The  comments  on  the 
Steinberg- Hunter  review  reflect  the  general  opinion  that  computer  generated  designs  present  a  very 
promising  area  for  the  future.  This  is  certainly  true  in  regard  to  the  field  of  RSM  designs.  The 
'danger'  that  is  projected  by  concerned  academicians  is  even  more  acute  in  the  RSM  area.  Many 
(though  not  all)  methods  of  generating  experimental  designs  employ  special  criteria  that  do  not 
address  robustness  notions  such  as  model  misspecification,  anticipation  of  outliers,  errors  in  control, 
extrapolation,  etc.  Model  misspecification  is  the  most  crucial  of  these  since  one  often  presumes 
from  the  outset  that  the  model  is  empirical  and  is,  indeed,  an  approximation.  All  of  this  suggests 
that  preparing  for  the  nonideal  and  very  detailed  interactive  properties  are  high  priorities  in  a  com¬ 
puter  generated  RSM  design.  We  would  like  to  underscoie  the  opinion  that  construction  of  the 
'export  system'  in  computer  driven  design  must  be  a  product  of  considerable  deliberation  and  pa¬ 
tience. 

There  appears  to  be  some  need  for  the  development  of  nonparametric  techniques  in  RSM. 
Most  of  our  analytic  procedures  depend  on  a  model.  The  use  of  model-free  techniques  would  avoid 
the  assumption  of  model  adequacy  or  low  order  polynomial  approximations  and,  in  particular,  the 
imposed  symmetry  associated  with  a  second  degree  polynomial.  The  use  of  nonparametric  kernel 
estimation  shows  some  promise. 

Analytic  methods  need  to  be  developed  for  the  case  of  repeated  measures  and  other  situations 
where  correlated  observations  are  encountered.  This  is  a  frequent  occurrence  in  biological  appli¬ 
cations. 
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